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ABSTRACT 

The  Naval  Weapons  Support  Center  is  planning  to  implement  a  bonus  system  to 
improve  the  reliability  of  pyrotechnic  devices.  The  measure  of  effectiveness  that  they 
wish  to  use  to  determine  how  to  award  bonuses  is  the  reliability  of  pyrotechnic  devices. 
The  data  available  to  estimate  this  reliability  is  based  on  the  current  sampling  inspection 
plan  in  which  devices  are  tested  in  different  environments.  The  models  which  include 
both  dependence  and  independence  assumptions  between  the  outcomes  of  these  tests  are 
implemented  and  estimates  of  overall  reliability  along  with  95  %  lower  confidence  bound 
are  obtained.  The  95  %  lower  confidence  bounds  are  found  by  bootstrapping.  Using 
these  estimates,  models  for  making  the  decision  to  award  bonuses  are  discussed  and 
studied  using  Monte  Carlo  simulation  . 
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I.     INTRODUCTION 

A.     BACKGROUND 

Nations  spend  a  lot  of  money  to  establish  a  strong  defense  network.  It  is  essential 
that  nations  buy  reliable  weapons  and  ammunition  from  the  contractors.  To  ensure  re- 
liability, contracts  must  include  a  lot  acceptance  sampling  plan  that  specifies  the  minimal 
acceptable  quality. 

The  Naval  Weapons  Support  Center  purchases  pyrotechnic  devices.  Unless  other- 
wise specified  in  the  contract,  the  supplier  is  responsible  to  see  that  his  devices  meet  all 
inspection  requirements  as  specified.  The  inspection  requirements  are  particular  to 
characteristics  of  each  type  of  device  and  are  specified  on  the  reference  drawings  and 
supplemental  quality  assurance  provisions  of  the  contract.  Any  testing  that  needs  to  be 
done  on  these  devices  is  explained  in  this  contract. 

When  nations  buy  weapons  and  ammunition  from  the  same  contractor,  they  would 
like  the  quality  to  improve  over  time.  As  contracts  are  now  written,  contractors  need 
only  to  satisfy  the  requirements  of  the  sampling  inspection  plan  for  lot  acceptance. 
Under  such  contracts,  contractors  have  no  incentive  to  improve  the  quality  of  items  they 
provide.  For  this  reason,  to  improve  quality,  The  Naval  Weapons  Support  Center  has 
decided  to  implement  a  bonus  system.  The  contractor  will  be  awarded  a  bonus  if  the 
result  of  the  sampling  inspection  exceeds  the  minimum  requirements  for  lot  acceptance. 

The  Naval  Weapons  Support  Center  will  begin  to  implement  a  bonus  system  for 
pyrotechnic  devices  in  FY91.  The  data  available  to  make  the  decision  whether  to  award 
a  bonus  is  based  on  the  current  sampling  inspection  plan.  This  plan  is  a  series  of  de- 
structive tests  in  different  environments.    The  purpose  of  this  thesis  is  to  provide  the 
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Naval  Weapons  Support  Center  with  guidance  for  implementing  a  bonus  system  for 
pyrotechnic  devices.  Once  implemented,  the  pyrotechnic  bonus  system  will  serve  as  a 
prototype  for  bonus  systems  for  other  devices. 

B.     PYROTECHNIC  DEVICE  RELIABILITY 

A  pyrotechnic  device  is  a  chemical  and  grenade  ammunition.  There  are  three  cate- 
gories of  pyrotechnic  devices.  The  three  categories  are:  Aerial  display,  Surface  display 
and  Grenades  fRef.  1:  pp.  1-2].  Samples  of  the  Pyrotechnic  devices  are  exposed  to 
various  environments  and  then  activated.  The  criteria  for  successful  activation  depends 
on  the  type  of  device  : 

1.  Aerial  display  (Ground  signals,  flares,  airburst  simulators,  sub  signals,  signal  kits, 
etc.) 

"Successful  activation  means  that  the  item  will,  after  simulating  user  environment, 
successfully " 

•  launch, 

•  have  proper  separation  /  signal  ignition, 

•  reach  desired  altitude  at  correct  angle, 

•  have  proper  parachute  deployment, 

•  have  proper  display  color, 

•  have  proper  display  time, 

•  have  no  subsequent  interference  with  next  item. 

2.  Surface  (Ground  or  Water)  display  (flare,  hand-held  signals,  smoke  /  illume 
grenades,  simulators,  etc.) 

"Successful  activation  means  that  the  item  will,  after  simulating  user  environment, 
successfully  ...." 


•  have  signal  ignition  with  proper  display, 

•  have  proper  display, 


•  have  proper  display  time, 

3.  Grenades  (fragmentation,  defensive,  white  phosphors,  etc.) 

Successful  activation  means  that  the  item  will,  after  simulating  user  environment, 
successfully  ...." 

•  function  (high  order  detonation  following  delay), 

•  have  proper  dissipation  of  payload, 

•  have  completely  consumed  payload. 

C.     PYROTECHNIC  DEVICE  RELIABILITY  PROBLEM 

Samples  from  any  large  lot  of  pyrotechnic  devices  submitted  by  a  manufacturer  must 
activate  after  exposure  to  different  environments.   These  are; 

1.  Manufacturer  Environment, 

2.  Temperature  and  Humidity  Environment, 

3.  Vibration  Environment, 

4.  Altitude  Environment. 

All  items  tested  are  subjected  to  the  manufacturer  environment.  However  items  are  only 
subjected  to  one  of  the  three  remaining  environments  :  Temperature  and  Humidity, 
Vibration  or  Altitude. 

The  sampling  plan  consists  of  using  four  distinct  samples  from  a  lot  that  can  be  as- 
sumed (approximately)  statistically  independent.  The  items  tested  in  each  sample  are 
also  assumed  to  be  independent.   According  to  the  sampling  plan; 

•  20  items  are  subjected  to  the  Manufacturer  Test, 


•  20  items  are  subjected  to  both  the  Temperature  and  Humidity  Test  and  Man- 
ufacturer Test, 

•  32  items  are  subjected  to  both  the  Vibration  Test  and  Manufacturer  Test, 

•  20  items  are  subjected  to  both  Altitude  Test  and  Manufacturer  Test. 

A  total  of  92  items  are  tested. 

Acceptance  criteria  for  each  test  are  : 

1.  Manufacturer  Test  :  Of  the  20  items;  if  no  more  than  1  fails  to  activate,  the  lot 
passes. 

2.  Joint  Temperature  and  Humidity  and  Manufacturer  Test  :  Of  the  20  items;  if  no 
more  than  1  fails  to  activate,  the  lot  passes. 

3.  Joint  Vibration  and  Manufacturer  Test  :  Of  the  32  items;  if  no  more  than  2  fails 
to  activate,  the  lot  passes. 

4.  Joint  Altitude  and  Manufacturer  Test  :  Of  the  20  items;  if  no  more  than  1  fails  to 
activate,  the  lot  passes. 

The  number  of  failures  for  these  tests  will  be  summarized  by  the  vector  : 

(   FOM,   FOTH,   FOV,   FOA    )  (1.2) 

where 

•  FOM  represents  the  number  of  failures  after  the  manufacturer  test, 

•  FOTH  represents  the  number  of  failures  after  both  the  temperature  and  humidity 
and  manufacturer  tests, 

•  FOV  represents  the  number  of  failures  after  both  the  vibration  and    manufacturer 
tests, 

•  FOA  represents  the  number  of  failures  after  the  altitude  and  manufacturer  tests. 


For  example, 

(   I,    1,   2,    1    ).  (1.2) 

represents,  1  failure  in  manufacturer  test,  1  failure  in  the  joint  temperature  and  humidity 
and  manufacturer  test,  2  failures  in  the  joint  vibration  and  manufacturer  test,  and  1 
failure  in  the  joint  altitude  and  manufacturer  test.    It  is  also  the  maximum  number  of 
failures  in  each  test  that  still  leads  to  lot  acceptance. 

The  marginal  distributions  of  the  number  of  failures  in  each  of  the  above  tests  are 
modeled  by  binomial  distributions.  There  are  24  possible  realizations  of  the  sampling 
inspection;  ranging  from  the  best  case  with  no  failures  to  activate,  to  the  worst  case  with 
1,  1,  2  and  1  devices  failing  to  activate  in  these  tests  manufacturer,  temperature  and 
humidity,  vibration  and  altitude  respectively.   These  cases  are  tabulated  in  Table  1 . 

To  award  bonuses  we  need  one  measure  of  effectiveness  for  pyrotechnic  devices  that 
can  be  estimated  from  the  available  data.  Ideally,  this  measure  is  the  reliability  of  the 
device.  However,  because  each  of  the  tests  are  destructive,  there  is  no  one  natural  defi- 
nition of  reliability  for  these  devices.  To  be  on  the  conservative  side,  we  define  the  reli- 
ability of  a  device  to  be  the  probability  that  the  device  will  activate  after  exposure  to  all 
of  the  environments.  It  is  not  easy  to  estimate  this  reliability  from  the  sampling  plan 
data.  This  data  is  incomplete  in  the  sense  that  we  have  limited  information  about  the 
joint  probability  of  activation  after  exposure  to  more  than  one  environment.  To  try  to 
compensate  for  this  lack  in  the  data,  we  will  use  models  for  the  joint  distribution  of 
(  FOM,  FOTH,  FOV,  FOA  )  that  specify  particular  types  of  dependence  between  the 
events  that  a  device  activates  after  exposure  to  different  environments. 


Table  1.     POSSIBLE  CASES 


CASE  #'S 

FOM 

FOTH 

FOV 

FOA 

1 

0 

0 

0 

0 

2 

0 

0 

0 

1 

3 

0 

0 

1 

0 

4 

0 

1 

0 

0 

5 

1 

0 

0 

0 

6 

0 

0 

1 

1 

7 

0 

1 

0 

1 

8 

0 

1 

1 

0 

9 

1 

0 

0 

1 

10 

1 

0 

1 

0 

11 

1 

1 

0 

0 

12 

0 

0 

2 

0 

13 

0 

1 

1 

1 

14 

1 

0 

1 

1 

15 

1 

1 

0 

1 

16 

1 

1 

1 

0 

17 

0 

0 

2 

1 

18 

0 

1 

2 

0 

19 

1 

0 

2 

0 

20 

1 

1 

1 

1 

21 

0 

1 

2 

1 

22 

1 

0 

2 

1 

23 

1 

1 

2 

0 

24 

1 

1 

2 

1 

Using  these  models  and  based  on  sampling  plan  data,  estimates  of  the  overall  reli- 
ability along  with  lower  confidence  bounds  are  obtained.  These  will  be  used  to  imple- 
ment the  bonus  system  for  pyrotechnic  devices.  We  compute  the  maximum  likelihood 
estimator  (  MLE  )  of  the  reliability  by  maximizing  the  appropriate  likelihood;  lower 
confidence  bounds  (  LCB  )  are  found  by  bootstrapping.  The  MLE's  are  computed  under 
both  independence  and  dependence  assumptions.  The  estimation  procedures  assuming 
independence  are  described  in  Chapter  II.  In  Chapter  III  we  incorporate  dependence 
by  fitting  a  Log  Linear  Model  to  our  data.  The  MLE's  from  Chapter  II  and  Chapter 
III  lead  to  inappropriate  results  for  this  for  this  problem;  thus  in  Chapter  IV  we  consider 
alternate  and  very  conservative  estimates  of  reliability.  Using  the  estimates  of  Chapter 
IV,  we  investigate  a  sequential  scheme  for  making  the  decision  to  award  bonuses  in 
Chapter  V.  The  results  of  simulations  are  presented  in  Chapter  VI.  Finally,  conclusions 
and  recommendations  are  given  in  Chapter  VII. 


II.     THE  ESTIMATION  OF  THE  MAXIMUM  LIKELIHOOD  ESTIMATOR 
(MLE)  OF  THE  RELIABILITY  WITH  INDEPENDENCE  ASSUMPTION 

A.     DEFINITIONS 

We  will  say  that  a  device  survives  environment  E  ,  if  it  is  still  potentially  capable  of 
activation  after  exposure  to  environment  E  .    Let, 

•  Ex  be  the  device  activates  after  exposure  to  Manufacturer  environment, 

•  £2  be  the  device  activates  after  exposure  to  Temperature  and  Humidity  environ- 
ment, 

•  £3  be  the  device  activates  after  exposure  to  Vibration  environment, 

•  £4  be  the  device  activates  after  exposure  to  Altitude  environment. 
We  define  the  reliability  of  device  as  below, 

r  =  p{e}  n  £2  n  £3  n  £*)•  (2.1) 

In  this  formula,  R  means  the  probability  that  a  device  activates  after  exposure  to  four 
environments.  We  will  estimate  R  for  each  of  the  24  cases  which  lead  to  lot  acceptance. 
Let  Q,  =  P(  El )  be  the  probability  that  device  activates  after  exposure  to  environ- 
ment i,  for  i  =  1,  2,  3,  4.  In  the  acceptance  sampling  plan  several  of  the  items  must 
activate  after  exposure  to  a  joint  manufacturer  and  another  environment.  To  avoid 
confusion  we  will  denote  tests  1  through  4  as  the  manufacturer  test,  the  joint  temper- 
ature humidity  and  manufacturer  test,  the  joint  vibration  and  manufacturer  test  and 
joint  altitude  and  manufacturer  test  respectively. 


Let 

*i  -  P(Ei)-  Qi  (2-2) 

and  let 

*,  =  P(E,n  %)  .  (2.3) 

Here  /?,  is  the  probability  that  a  device  survives  test  i  for  i  =  1,  2,  3,  4.  The  simplest 
model  is  to  assume  that  £,,  E^,  £,,  £4  are  independent.  If  we  assume  that  Elt  .  .  .  ,  £4 
are  independent  then 

Ri  =  Qx  Qt  (2-4) 

for  i  =  2,3,4  and  the  reliability  of  device  is, 

R=  QxQi  Qi  Qa  ■  (2-5) 

B.     THE  LIKELIHOOD  EQUATION 

Let 

•  X,  be  the  number  of  devices  that  activate  after  test  i, 

•  nt  be  the  number  of  items  given  test  i. 

Then  X,  is  binomial  with  parameters  R,  and  n,  for  i=  1,  2,  3,  4.  Under  the  assumption 
of  independence  the  joint  likelihood  function  of  observing  Xx  ==  jcx  ,  . . . ,  XA  =  x4  is 


L(xl,x2,x2,x4\Rl,R2,R3,R4)  =         P  J  R?  (  1  -  Rj-  '*■  (2.6) 


with  constraints  : 

0   <    #!    <    1, 

0  <  R2  <  Ru 

0  <  R2  <  /?,, 

0  <  R*  <  Rv 

Our  aim  is  to  maximize  this  likelihood  function  subject  to  the  constraints  that 
(Rlf  R2,  R3,  R4)eS  where 

5  =  { (  Rx  R2  R3  #4 ) :  0  <  Rx  <   1,    0  <  Rt  <  Rx      i  =  2,  3,  4  }  . 
From  the  equation  (2.6),  we  see  that  maximizing  L  is  equivalent  to  maximizing 

4 

/  =  £{(*,  In*,)  +  (/i,  -  xt)  ln(l-/?,)},  (2.7) 

i=l 

where  the  constant  multipliers  I      J    for  i  =  1,  2,  3,  4  have  been  dropped  (  because  they 
do  not  effect  the  maximization  procedure  )  and  the  natural  logarithm  of  L  is  taken. 

We  first  show  that  /  is  a  concave  function.  To  show  that  /  is  a  concave  function, 
we  can  show  -  /  is  a  convex  function.  According  to  Theorem  3.3.6  [Ref.  2:  p.  92], 
by  looking  at  its  Hessian  matrix,  we  can  learn  whether  function  is  convex  or  not.  If  its 
Hessian  matrix  is  positive  semi-definite  at  each  point  S  then  function  /  is  convex.  To 
create  the  Hessian  matrix,  we  must  calculate  partial  derivatives  of  the  function   -   / , 


ft  I  Xi  n,    —     X: 

+  t — ir  (2-8) 


dRi  Ri  1  -  Rj 
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r>2  /  xi  ni  ~   xi 

a^;2        /??        i  -  /?,2 


-,2  . 

=  0  i*j.  (2.10) 


a  /?,  Rj 


Then,  the  determinant  of  the  Hessian  is: 


njc.-  n,  —  x> 


<=i 


1   ■  •   R? 


Clearly  we  can  see  that  for  0  <  R,  <   1  i  =  1,  2,  3,  4,    |  H  |  is  always  positive.  Because 
-    /   is  continuous,  this  implies  that   -   /  is  a  convex  function  on  S.   As  a  result  of  this, 
/  is  a  concave  function. 

We  note  that  with  the  constraints  on  the  probabilities  /?, ,  .  .  .  ,  /^  ,  there  does  not 
in  general  exist  a  closed  form  solution  to  MLE.  However,  with  only  24  realizations  of 
j:,  ,  . . .  ,  x4  of  interest,  the  estimated  reliabilities  for  these  24  cases  can  be  found  with 
some  rather  tedious  but  straight-forward  computations. 

C.     COMPUTING  THE  MAXIMUM  LIKELIHOOD  ESTIMATORS  AND  LOWER 
CONFIDENCE  BOUNDS 

Because  /  is  concave  over  the  convex  set  S,  if  the  maximum  occurs  in  the  interior 
of  S,  it  is  a  unique  maximum  and  is  given  by 


Rt  =  it  (2-12) 


for  i  =  1,2,  3,  4.  where 
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A  A 


0  <  /?!  <  1     ,   0  <  Ri  <  #, 

for  i  =  2,  3,  4. 

In  this  case  the  MLE's  for  Q,  i  =   1,  2,  3,  4  are 


<2i  =  -^  (2.13) 


A  A7: 

Qi  =   1FT  (2-14) 


A7 


for  i  =  2,  3,  4.    Finally  we  can  estimate  the  reliability  of  the  device  as 

R  =  QxQiQiQa-  (2-15) 


Table  2  summarizes  the  cases  for  which  the  MLE's  can  be  found  using  (2.12)  - 
(2.15).  When  the  maximum  falls  on  the  boundary  of  S,  there  is  no  explicit  expression 
for  the  MLE  of  (/?,,  . . . ,  R+ ) .  These  are  the  cases  with  the  exception  of  (  1  12  1) 
which  have  one  failed  item  in  the  manufacturer  test.  This  implies  that  /?,  will  be  less 
than  1.0  .  To  find  the  MLE,  we  find  the  maximum  of  (  Ry ,  . . ,  R+)  on  each  of  the 
boundaries,  compute  /  for  each  of  these  and  let  the  MLE  be  the  one  with  the  largest 
value  of  /  .  It  is  clear  that  in  most  cases  several  of  the  boundaries  can  be  eliminated 
from  consideration,  simplifying  computation  considerably. 
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Table  2.     EASY  AND  HARD  CASES 

EASY  CASES 

HARD  CASES 

(0        0        0        0     ) 

(1        0        0        0     ) 

(0001) 

(10        0        1) 

(0        0        1        0     ) 

(1        0        1        0     ) 

(0        1        0        0     ) 

(1         10        0     ) 

(001         1     ) 

(1        0        1         1     ) 

(0        1        0        1     ) 

(1         1        0        1     ) 

(0        1        1        0     ) 

(1         1         1        0     ) 

(0        0        2        0     ) 

(10        2        0) 

(0        1        1        1     ) 

(1      1      1      1    ) 

(0        0       2        1     ) 

(10        2        1) 

(0        1        2        0     ) 

(1         12        0) 

(0        1        2        1     ) 

(112        1) 

D.     EXAMPLES 
1.     Easy  Case 

In  this  example  our  failure  vector  is, 


(0121). 
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According  to  this  failure  vector,  we  can  write  our  likelihood  equation  by  using  equation 
(2.6)  to  get 

L  =  R?{  1  -  /?,  )°/?219(  1  -  R2  )'  R?{  1  -  R2?R\\\  -  R<  )'     (2.16) 

0     <      /?!      <      1      ,       0     <£,.     </?!     , 

for  i  =  2,  3,  4.   From  the  likelihood  above; 


A  *i  20 

R<  -  -sr  -  20  -  1000° 


A  X->  1Q 

^  -  T5"  -  lr  -  °-9500 


*3  =   "*    =     32    =  °9375 


^  =    "S"   =    "20"    =  °-9500'  (2>17) 


Results  imply  that  all  R,'s  for  i  =  2,  3,  4  are  between  0.0  and  Ru  which  means  that 
constraints  are  met.   Now  we  can  estimate  Q,'s  using  equations  (2.13)  and  (2.14). 

Qx  =  1.0000 

Q2  =  °-9500  =  0.9500 

q2   =  JL22Z1  _  0.9375  (218) 

&  =  °-9500  =  0.9500 . 

And  finally,  reliability  of  the  device  can  be  estimated  by  using  equation  (2.15) 
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R  =  (1.0000)  (0.9500)  (0.9375)  (0.9500)  =  0.84609375.  (2.19) 

2.     Hard  Case 

In  this  example  our  failure  vector  is 

(10        2        1) 

According  to  this  failure  vector,  we  can  write  our  likelihood  equation  by  using  equation 
(2.6)  as 

L  =  R\9  (  1  -  Rx  )}  R?  (  1  -  R2f  R?  (  I  -  R3)2  Rl9(  I  -  R4  )'     (2.20) 
0    <    /?,    <    1    ,    0    <Rt   <RX 
for  i  =  2,  3,  4.    From  the  likelihood  above; 


— -   =   —   =  0  9500 
«i  20 


xi  20 


*3  30 

—   =   —   =  0  9375 
«3  32 


i   =   IS.   =  0.9500.  (2.21) 


JC  JC 

As  you  see  from  the  above  -=-   >   -=-  thus,  the  MLE  does  not  lie  in  the  interior  of  S. 

"2  "1 

We  begin  by  considering  the  boundary 
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/?,  —  R2  —  R3  —  /?4 . 


The  likelihood  equation  on  this  boundary  is, 


L  =  /?f8  (  1  -  Rx  f  (2.22) 


0    <    Rx    <    1 


From  the  likelihood  equation, 


A  A  A  A  CO 

Rx  =  R2  =  R3  =  R4  =    J5-    =  0.9565 


Then  value  of  the  likelihood  with  these  estimated  /?,'s  is 


L  =  (  "I"  )88  (  ~k  )4  =  °-714  X   10~7  "  (2'23) 


For  the  boundary 

Rx  =  R2  =  /?3  ,  /?4  <  R\  , 
the  likelihood  equation  is, 


Z,  =  tff9  (  1  -  K,  )3  Rl9  (  1  -  ^  y  (2.24) 


0    <    /?,    <    1  ,     0</?4<^i. 
From  the  likelihood  equation, 


A  A  A  (LQ  A  1Q 

R^  =  R2=  R2  =   ^.   =  0.9583  R4=   -^.   =  0.9500 , 
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and  the  value  of  the  likelihood  with  these  estimated  R's  is 


69     ,69,      3      ,3,     19     ,19,       1 


L  =  (    72    )M(    12    )3(    20    )19(  "20    }'  =  °-723  X   10"7,  (2-25) 


For  the  boundary 

/?,  =  R2  =  R4  ,  R-2  <  R\ 
the  likelihood  equation  is, 


L  =  R\*  (  1  -  Rx  f  R?  (  1  -  R3  f  (2.26) 


0    <    Rx    <    1        0    <    R3    <  R} 
From  the  likelihood  equation, 


A  A  A  <C 

fl,  =  /?2  =  /?4  =    -gL    =  0.9666 


/?3  =   -£-   =  0.9375 . 
J         32 


Then  value  of  the  likelihood  with  these  estimated  R,'s  is 


L  =  (  "60"  )58(  ~k  )2(  ~32~  )3°(  ~32~  )2  =  °-876  X   10~7  •  (2'27) 


For  the  boundary 

/?,  =  R3  =  R4  ,  R2  <  R\  , 
the  likelihood  equation  is, 
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L  =  R™  (  1  -  /?,  )4  R?  (2.28) 


0    <    Rx    <    1       0    <    /?2    <  /?,  . 

Clearly  L  is  maximized  on  the  boundary  of  the  constraints  (2.28),  i.e.  on  another  of  the 
boundaries  of  S  thus  we  can  eleminate  this  case  from  consideration.    For  the  boundary 

Rx  =  R2  ,  R3  <  /?,  ,  R4  <  Rl  , 
the  likelihood  equation  is, 

L  =  /?r  (  1  -  /?,  )l  R™  (  1  -  R2  )2  Rl9  (  1  -  #4  )'  (2-29) 

0    <    Rt    <    1  ,  0    <    R2    <  #,  ,  0    <    R4    <  /?!  . 
From  the  likelihood  equation, 


A  39  30  A  19 

R    =  R    =   J2-  =  0.9750       R3  =   ^-   =  0.9375        R^  =   -£r   =  0.9500 
z  40  32  20 


Then  value  of  the  likelihood  with  these  estimated  R's  is 


L  -  <  lr  A  io )l(  1" )30(  ■£ )2(  lr )l9(  ■»  >'  =  °"° x  10"7<2-3°) 


Quick  inspection  reveals  that  the  remaining  boundaries  can  be  eliminated  from  consid- 
eration. 

After  boundary  analysis,  we  can  see  that  equation  (2.30)  gives  us  the  maximum 
likelihood  value.   Thus  the  MLE's  for  /?, ,  R2 ,  R3 ,  /^  are 


/?,  =  0.9750 
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#2  =  0.9750 


R,  =  0.9375 


RA  =  0.9500 


Now  we  can  estimate  <2,'s  with  equations  (2.13)  and  (2.14) 


<2,  =  0.9750 

&  =  Wm  =  L0000 


„  0.9375         0Q,  s  (2.31) 

&  =  O9750"  ~  °-%15 

_  0.9500 

^4  "  0.9750         u-y/4^- 


And  finally  reliability  of  the  device  can  be  estimated  by  using  equation  (2.15) 

R  =  ( 0.9750 )  (  1.0000 )  ( 0.9615  )  ( 0.9744  )  =  0.91346100  .  (2.32) 

We  computed  MLE's  with  a  FORTRAN  program  MLEA  given  in  Appendix 
A.  After  24  replications  of  this  program,  we  can  estimate  R,  's  for  i  =  1,  2,  3,  4  for  each 
of  the  24  outcomes  of  the  sampling  plan  that  lead  to  lot  acceptance.  These  are  given  in 
Table  3,  R  is  given  in  Table  4. 

We  compute  lower  confidence  bounds  for  each  possible  case  by  bootstrapping. 
The  bootstrap  can  be  used  to  produce  approximate  confidence  intervals  in  an  automatic 
way.  There  are  several  ways  to  set  approximate  confidence  intervals  with  bootstrapping. 
These  are  the  percentile  method,  the  standard  method,  the  bias-corrected  percentile 
method  and  the  nonparametric  method  [Ref.  3  :  pp.  67-70].  Let  6  be  an  unknown  pa- 
rameter with  estimator  6.    To  bootstrap,  samples  are  generated  using  6  in  place  of  the 
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unknown  parameter  6.  For  each  of  these  samples  an  estimator  8  is  computed.  We 
define  G(s)  to  be  the  parametric  bootstrap  cumulative  distribution  function  of  6',  i.e. 
G  ( 5 )  is  the  empirical  distribution  of  the  8"s.  All  methods  mentioned  above  use 
percentiles  of  G  to  define  the  confidence  interval.  They  differ  in  which  percentiles  are 
used.   The  percentile  method  was  used  in  our  calculations. 

If  we  use  the  notation  6  ( a  )  for  the  level  a,  endpoint  of  one  sided  lower  confi- 
dence interval  for  9  ,  then 

/>[0(a)<0]  =  a.  (2.33) 

An  estimate  of  8  (  a  )  from  the  bootstrap  cumulative  distribution  is  given  by 


0(a)  m  G~]  (a).  (2.34) 


To  get  LCB's  for  R,  for  each  of  the  24  realizations  of  the  sampling  plan  ,  we 
generate  bootstrap  samples  of  random  failure  vectors,  in  which  failures  come  from  in- 
dependent binomial  distribution  with  parameters  (  n, ,  R, )  for  i  =  1,  2,  3,  4.  This  is  done 
using  the  FORTRAN  program  RANVEC  [Ref.  4)  in  Appendix  B.  We  generate  5000 
failure  vectors  for  each  case.  And  then  we  estimate  R  s  from  each  of  the  5000  failure 
vectors  by  the  means  of  the  program  MLEA.    The  next  step  is  to  compute  the  order 

A  A  A 

statistics  of  R  's  from  /?,  to  Rsooo.  We  get  the  parametric  bootstrap  cumulative  distrib- 
ution function  (  i.e.  the  emprical  distribution  of  /J, ,  . . . ,  R^on  )  with  this  computation. 
This  is  done  using  the  FORTRAN  program  SORT  in  Appendix  C. 

Finally  we  can  get  the  95  %  lower  confidence  bound  using  equation  (2.34)  from 
this  routine.  Reliabilities  and  95  %  lower  confidence  bounds  are  listed  in  Table  4.  Re- 
sults in  Table  4  are  given  in  descending  order  and  R's  and  LCB's  are  not  ordered  as  we 
expected  them  to  be.    For  example,  the  failure  vector  that  has  the  maximum  number  of 
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failure  in  each  test  is  in  the  middle  of  the  table  with  respect  to  R.  It  has  also  a  bigger 
95  %  LCB  than  the  failure  vector  that  has  a  total  of  3  failed  items  in  each  of  the  joint 
tests.  These  results  are  counter-intuitive  because  we  expect  that  more  failures  indicate 
a  lower  overall  reliability.  This  is  reasonable  because  it  is  likely  that  a  device  that  is 
poorly  constructed  is  more  likely  to  fail  in  any  of  the  four  environments.  From  the  re- 
sults in  Table  4,  it  is  clear  that  an  attempt  to  model  dependence  must  be  made  in  order 
to  get  believable  estimates  of  R. 
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Table  3.  TEST  PROBABILITIES  (  A  ) 

CASE 

R> 

R> 

R3 

K 

(00  00) 

1.0000000 

1.0000000 

1.0000000 

1.0000000 

(000  1) 

1.0000000 

1.0000000 

1.0000000 

0.9500000 

(00  10) 

1.0000000 

1.0000000 

0.9687500 

0.1.00000 

(0  100) 

1.0000000 

0.9500000 

1.0000000 

1.0000000 

(10  00) 

0.9891304 

0.9891304 

0.9891304 

0.9891304 

(0011) 

1.0000000 

1.0000000 

0.9687500 

0.9500000 

(0101) 

1.0000000 

0.9500000 

1.0000000 

0.9500000 

(0110) 

1.0000000 

0.9500000 

1.9687500 

1.0000000 

(1001) 

0.9861111 

0.9861111 

0.9861111 

0.9500000 

(10  10) 

0.9833333 

0.9833333 

0.9687500 

0.9833333 

(1  100) 

0.9861111 

0.9500000 

0.9861111 

0.9861111 

(0020) 

1.0000000 

1.0000000 

0.9375000 

1.0000000 

(0111) 

1.0000000 

0.9500000 

0.9687500 

0.9500000 

(10  11) 

0.9750000 

0.9750000 

0.9687500 

0.9500000 

(110  1) 

0.9807692 

0.9500000 

0.9807692 

0.9500000 

(1110) 

0.9750000 

0.9500000 

0.9687500 

0.9750000 

(002  1) 

1.0000000 

1.0000000 

0.9375000 

0.9500000 

(0  120) 

1.0000000 

0.9500000 

0.9375000 

1.0000000 

(1020) 

0.9833333 

0.9833333 

0.9375000 

0.9833333 

(1111) 

0.9615384 

0.9500000 

0.9615384 

0.9500000 

(0121) 

1.0000000 

0.9500000 

0.9375000 

0.9500000 

(1021) 

0.9750000 

0.9750000 

0.9375000 

0.9500000 

(1 120) 

0.9750000 

0.9500000 

0.9375000 

0.9750000 

(112  1) 

0.9500000 

0.9500000 

0.9375000 

0.9500000 
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Table  4.     RELIABILITIES  AND  95  %  LOWER  CONFIDENCE  BOUNDS 

(A) 


FAILURE 
VECTOR 

A 

R 

FAILURE 
VECTOR 

95  %  LCB 

(0  0  00) 

1.0000000 

(000  0) 

1.0000000 

(1000) 

0.9891304 

(00  10) 

0.9062500 

(00  10) 

0.9687500 

(1000) 

0.9000000 

(10  10) 

0.9687499 

(0020) 

0.8750000 

(000  1) 

0.9500000 

(10  10) 

0.8550000 

(0  100) 

0.9500000 

(000  1) 

0.8500000 

(1001) 

0.9499999 

(0  100) 

0.8500000 

(1100) 

0.9499999 

(1001) 

0.8282812 

(10  11) 

0.9439102 

(1100) 

0.8258822 

(1110) 

0.9439102 

(0110) 

0.8234775 

(1111) 

0.9385999 

(00  1  1  ) 

0.8234375 

(0020) 

0.9375000 

(1020) 

0.8125000 

(1020) 

0.9374999 

(10  11) 

0.8015624 

(112  1) 

0.937499 

(0101) 

0.8000000 

(00  11) 

0.9203125 

(1110) 

0.8000000 

(0  110) 

0.9203125 

(0  120) 

0.7875000 

(110  1) 

0.9201959 

(002  1) 

0.7749999 

(1021) 

0.9134614 

(110  1) 

0.7749018 

(1120) 

0.9134614 

(1111) 

0.7649999 

(0101) 

0.9025000 

(1021) 

0.7647058 

(002  1) 

0.8906249 

(1120) 

0.7614843 

(0  120) 

0.8906249 

(112  1) 

0.7505192 

(0111) 

0.8742968 

(0111) 

0.7499999 

(0  12  1) 

0.8460937 

(0  12  1) 

0.7171874 
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III.     LOG  LINEAR  MODEL  WITH  DEPENDENCE  ASSUMPTION 

A.     BACKGROUND 

Modeling  the  outcomes  of  the  various  environmental  tests  as  independent  is  clearly 
inappropriate.  However,  the  nature  of  the  acceptance  sampling  plan  makes  it  impossible 
to  estimate  the  reliability  of  the  device  (  i.e.  the  probability  that  it  would  activate  after 
exposure  to  all  four  environments  )  without  some  assumptions  about  the  dependence 
between  outcomes  of  various  tests.  Thus,  our  second  approach  is  to  model  the 
pyrotechnic  device  reliability  with  a  log  linear  model. 

The  results  of  our  test  series  create  a(2x2x2x2)    contingency  table.    Let    p 
represent  passing  a  test  and   f    represent  failing  a  test. 

Table  5.     CONTINGENCY  TABLE  STRUCTURE   FOR  TEST  SERIES 


Passed 
Manufacturer  Test 

Failed 
Manufacturer  Test 

Passed 

Tem.&Hum. 

Test 

Failed 

Tem.&Hum. 

Test 

Passed 

Tem.&Hum. 

Test 

Failed 

Tem.&Hum. 

Test 

Passed 

Vibration 

Test 

Passed 

Altitude 

Test 

Mpppp 

MPfPP 

MfPPP 

MJTpp 

Failed 

Altitude 

Test 

MPPPf 

Mpfpf 

MfPPf 

Mffpf 

Failed 

Vibration 

Test 

Passed 

Altitude 

Test 

MPPfP 

MpJTp 

MfpfP 

MmP 

Failed 

Altitude 

Test 

MppIT 

MpJJT 

MfPff 

MM 
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The  results  of  our  tests  can  be  thought  of  as  censored  data  from  a  hypothetical 
(2X2X2X2)  contingency  table  (See  Table  1) .  The  frequency  in  each  cell  of  this  table 
is 

Mm  ,        ( i,  j,  k,  I)  e  {  p,  ff  , 

the  number  of  devices  out  of  92  which  would  have  result  i  in  environment  1  (manufac- 
turer), result  j  in  environment  2  (temperature  and  humidity  alone),  result  k  in  environ- 
ment 3  (vibration  alone),  result   /  in  environment  4  (altitude  alone). 

This  is  a  hypothetical  table,  because  if  a  device  was  exposed  to  all  four  environments 
and  then  failed  to  activate,  there  would  be  no  way  to  discern  which  combination  of  the 
four  environments  caused  failure.  The  data  from  the  acceptance  sampling  plan  can  be 
thought  of  as  censored  data  from  such  a  (2X2X2X2)  contingency  table.  As  an 
example,  a  device  that  is  given  just  the  manufacturer  test  belongs  in  one  of  the  cells 
(  p,  j,  k,  / )  where  ( j,  k,  / )  e  { p,  /}3,  because  it  is  not  clear  what  would  have  happened 
to  it  had  it  been  exposed  to  the  other  three  environments. 

Using  a  log  linear  model  under  independence,  the  expected  value  of  each  cell  fre- 
quency is 

where   i,  j,   k,     /    =    p,  f  and 

4  =  -  4  (3.2) 

for  i  =  1,2,  3,  4.  Here  the  parameter  \x  represents  the  overall  effect  and  the  parameters 
X'p  represent  the  effect  of  passing  environment  i.  It  is  simple  to  show  that  the  log  linear 
model  (3.1)  is  equivalent  to  the  independence  of  the  environments  [Ref.  5:  pp.  25-46]. 
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On  a  logarithmic  scale,  the  independence  relation  is  equivalent  to  the  additive  re- 
lationship.  As  an  example, 


log  Mpppp  =  n  +  Alp  +  X2p  +  4  +  ij .  (3.3) 


Log  linear  models,  which  take  into  account  dependence  include  extra  interaction  terms. 
Because  of  the  extreme  amount  of  censoring,  we  will  only  consider  models  with  two  way 
interaction  terms  of  the  type. 

■  i'  j  i2  ^  j3  ^  J4  _■_  j12  _■_  in  _■_  i14  j.  j23 _■_  i24  _i_  j34 
E  T  M     "1  =  e^      (      J      *  '  *   '*   "   ■**  ''   *'        (3.4) 

where   i,  j,  k,   /    =  p,  f, 

4  =    -  4  (3.5) 

for  i  =  1,  2,  3,  4  and 

4  -  4  -  -  4  -  -  4  <3-6) 

for  i,  j,  =  1,2,  3,  4.  As  an  example; 

i    i>    .    j2   ■    j3   .    J4  i    i'2   .     .I3   .    .14   .    ,23   .    ,24  j    ,34 

E\_M       ']  =  e't+'+'+'+'+''+»*»+»+»''.  (3.7) 

The  censored  data  does  not  tell  us  much  about  partial  association  between  pairs  of 
tests.  For  this  reason,  we  assume  that  partial  associations  between  any  pair  of  tests  are 
all  the  same.  Let  6  be  the  partial  association  between  tests.  Then  we  can  reformulate 
equation  (3.7)  as  follows: 

ELM^l-e"*'^''^^''^6.  (3.8) 
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The  MLE's  for  the  expected  number  of  devices  in  each  cell  (  or  equivalently  the 
MLE's  for  the  parameters  )  can  not  be  found  explicitly.  The  EM  Algorithm  (  Expecta- 
tion -  Maximization  )  will  be  used  to  approximate  the  MLE's  of  the  expected  number 
of  devices  in  each  cell. 

B.     EXPECTATION  MAXIMIZATION  (  EM  )  ALGORITHM 

A  general  method  of  maximum  likelihood  estimation  from  incomplete  data  is  the 
EM  Algorithm.  The  expectation  maximization  algorithm  is  an  iterative  procedure  where 
each  stage  consists  of: 

■  an  expectation  step  (  E  )  followed  by 

■  a  maximization  step  (  M  ). 

This  algorithm  is  generally  used  to  compute  maximum  likelihood  estimators  in  in- 
complete data  problems.    In  the  application  of  the  EM  algorithm  we: 

•  replace  missing  values  by  their  estimated  expected  values  given  the  incomplete  data 

•  estimate  parameters 

•  reestimate  the  missing  values  assuming  the  new  parameter  estimates  are  correct 

•  reestimate  parameters 

and  so  forth,  iterating  until  convergence  [Ref.  6:  pp.  127-141].   This  iterative  algorithm 
works  as  follows  in  pyrotechnic  device  problem. 
1.     Initialization 

The  algorithm  requires  initial  guesses  for  the  parameters  of  the  log  linear  model 
cell  frequency  in  contingency  table.  It  uses  initial  guesses  and  calculates  initial  cell 
probabilities.  We  will  get  our  initial  guesses  by  first  fitting  a  log  linear  model  under  in- 
dependence. 
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2.     Iterations 

■  Expectation  Step 

EM  Algorithm  estimates  expected  cell  frequencies  given  the  data  by  using  the 
most  current  estimates  of  the  parameters.  It  compares  these  estimated  expected  cell 
frequencies  with  the  previous  estimates.  If  the  differences  between  these  two  estimates 
are  small  enough,  then  the  EM  Algorithm  has  converged.  We  then  accept  the  final  es- 
timates of  cell  frequencies  as  MLE,  so  we  can  easily  calculate  cell  probabilities. 

■  Maximization  Step 

If  the  algorithm  has  not  converged,  then  it  starts  to  estimate  new  parameters  for 
the  log  linear  model  using  expected  cell  frequencies  from  the  E  step  as  if  there  were  ac- 
tual data  available.  Estimation  is  done  by  maximizing  the  likelihood  using  an  iterative 
Newton- Raphson  method.  The  estimated  parameters  from  this  step  are  than  used  in  the 
next  E  step  of  the  EM  algorithm. 

C.     CALCULATIONS 

We  apply  the  EM  algorithm  for  each  realization  of  the  failure  vector  (  FOM, 
FOTH,  FOV,  FOA  ).  During  the  maximization  step  of  the  EM  Algorithm,  we  use  the 
Newton-Raphson  procedure  which  is  described  by  SAS,  [Ref.  7:  pp.  190-212].  We  first 
describe  the  M  step  using  as  an  example  the  following  failure  vector. 

(112     0) 

■   Maximization  Step 

We  start  with  initial  guesses  (IG)  for  Mijkt's.  These  are  the  conditional  expected 
values  provided  by  the  previous  E  step.  For  example,  initial  guesses  for  the  failure  vector 
above,  are  shown  below. 
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IG 


M, 


pppp 


M 


pppf 


M 


ppfp 


M 


PPff 


M 


p/pp 


M 


Pfpf 


M 


pffp 


M 


Pfff 


If 


fppp 


M 


'fppf 


M 


'/P/P 


M 


fpff 


M 


f/pp 


M 


VfPf 


M, 


Wp 


M, 


W 


IG    = 


83.08 

0.1 
3.36 

0.1 
2.21 

0.1 

0.1 

0.1 
2.15 

0.1 

0.1 

0.1 

0.1 

0.1 

0.1 

0.1 


(3.9) 


From  these,  we  can  compute  the  proportion  of  observations  which  fall  into  each  cell  as: 


A 
P       = 


M, 


pppp 


92 


M, 


pppf 


92 


M, 


ppfp 


92 


M, 


ffff 


92 


0.9030 
0.0011 
0.0365 

0.0011 


(3.10) 


We  will  use  P  to  get  MLE's  for  the  parameters  of  the  log  linear  model.    Our  log  linear 
model  which  includes  two  way  interactions  can  be  written  as 


29 


log  Ppppp 
log  Ppppf 
log  Pppfp 

log  Pm 


1    1 
1    1 


1       1     6 
1     -1    0 


1    -1    -1    -1    -1    6 


y    -    p 


(3.11) 


where  yT  =  ( p.,  X\,  Xj,  X\,  X* )  is  the  parameter  vector  and  p  is  normalizing  constant 
required  by  the  restriction  that  the  probabilities  sum  to  1  and  where  P,Jkl  is  the  proba- 
bility corresponding  to  the  ( ij  k  /  )rt  cell  in  the  hypothetical  contingency  table.  We  can 
derive  from  equation  (3.1 1) 


log  Ppppp 

log  Ppppf 

l°ZPPP/P 

= 

logPM 

1       1     6 
1      -1    0 


-1    -1    -1    -1    6 


fi   -    3 


(3.12) 


where  fiT  =  ( Xlp,  X),  X],  X* )  is  the  parameter  vector  and  5  is  normalizing  constant  re- 
quired by  the  restriction  that  the  probabilities  sum  to  1.  Then,  to  use  the  SAS  proce- 
dure, we  rewrite  log  />„„,  log  Pm/t  .  .  .  ,  log  Pm  as  the  15  logits,  F,  =  ( log  PPJ  log  Pm), 
.  .  .  ,  Fls  =  ( log  Pnfp\  log  Pm)  so  that 


F(P)  = 


15 


10   0   0-1 
0    10   0-1 


0   0   0   0-1 


log 

p 

1  pppp 

log 

p 

1 pppf 

log 

p 

1 ppfp 

log 

PM 

(3.13) 
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where  P  =  (  Ppppp,  Pppp/,  . . . ,  Pm/).    Using  equations  above  the  following  result  is  ob- 
tained 


F(P) 


F, 


F 


15 


2    2    2    2     0 
2    2    2    0-6 

0    0    0    2-6 


(3.14) 


The  design  matrix  (X)  of  this  problem,  from  the  equation  (3.14)  is  as  follows: 


X  = 


2 

2 

2 

2 

0 

2 

2 

2 

0 

-6 

2 

2 

0 

2 

-6 

2 

2 

0 

0 

-8 

2 

0 

2 

2 

-6 

2 

0 

2 

0 

-8 

2 

0 

0 

2 

-8 

2 

0 

0 

0 

-6 

0 

2 

2 

2 

-6 

0 

2 

2 

0 

-8 

0 

2 

0 

2 

-8 

0 

2 

0 

0 

-6 

0 

0 

2 

2 

-8 

0 

0 

2 

0 

-6 

0   0   0    2-6 


(3.15) 


In  the  application  of  Newton-Raphson  Method,  we  use  the  variance  and  covariance 
matrix   S  of  F(  P ).    It's  inverse  is  given  by 


S'l{P)    = 


p      —Pi 

1  pppp     i  pppp 


P  V      P 

1  PPPP    *    l  PPPf 


P  v    P  P  —    P 

1  PPPP         l  PPPf      I  PPPf         l  PPPf 


PPPPP   X    PfffP     ~    PPPPf  X    Pfffp 


Ppppp   x    Pfffp 
Ppppf  X    PfffP 


PfffP   ~    PfffP 


(3.16) 
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The  first  estimates  of  parameters  for  the  Newton- Raphson  procedure  are  calculated 
as  follows: 

b0  =  [XTS~\{P)XY*  lXTS~]  (P)F{P)1-  (3-17) 

where  Sl  and  F  are  estimated  using  the  proportions  in  the  vector  P.  We  estimate  the 
reduced  logit  response  functions  using  the  equation 

FQ  =  Xb0.  (3.18) 

From  F0  we  can  compute  the  updated  estimates  of  Ppppp ,  Ppppf ,  . . . ,  Pm  Let 
n  =  (  n,  (  1  ),  n2  (  1  ),  .  . .  ,  ni6  (  1  ) )  be  the  vector,  which  contains  these  estimates  of 
the  cell  probabilities.   The  value  of  the  log  likelihood  evaluated  at  n  is 

16 

LHE  =  Yjxi  lQg  n(0  (319) 

where  x,  is  the  number  of  items  in  cell  i  with  respect  to  probability  in  n. 

We  estimate  parameters  iteratively  until  the  difference  between  last  estimate  and 
previous  is  small  enough.  At  each  iteration  we  update  the  inverse  of  variance  and 
covariance  matrix  with  probabilities  of  n  from  the  previous  iteration  and  we  do  fol- 
lowing matrix  computations. 


C  =  XTS~]{T\)X 

rr  *  -,  (3-20) 

G  =  Xr[92.0  x  (P  -  n)] 


Let  b,  be  the  next  estimate  of  parameters  /?  in  the  i  th  iteration.   Then  b,  as  follows: 


bt  =  bt_x  -  dC~x  G  (3.21) 
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where  5  <  1  is  a  constant  supplied  by  the  user. 

We  get  first  cell  expectations  from  initialization.  The  failure  vector,  which  is  defined 
at  page  5,  has  initial  guesses  for  cells  which  are  shown  in  equation  (3.9).  After  the  first 
application  of  Newton  -Raphson  Algorithm,  we  have  following  initial  cell  expectations 
(CE). 


CE   = 


82.4278259 

0.7278117 

3.1047230 

0.1208856 

2.0790358 

0.0809494 

0.3453168 

0.0592893 

2.0268116 

0.0789160 

0.3366428 

0.0578000 

0.2254283 

0.0387050 

0.1651089 

0.1250073 


(3.22) 


■   Expectation  Step 

The  conditional  expected  frequency  for  each  cell  expectations  are  calculated  using 
the  estimated  parameters  of  the  log  linear  model  from  the  previous  M  step.  Some  of  the 
cell  expectations  formulated  below 


E IMPPPP  |  (FOM,  FOTH,  FOV,  FOA)  ]  =  (20  - 


+  (20 
+  (32 
+  (20 


FOM)xP{pppp\p...) 
FOTH)xP(pppp\pp..) 
FOV)xP(pppp\p.p.) 
FOA)xP(pppp\p..p) 


(3.23) 


33 


ElMpfp/\(FOMfFOTH,FOV,FOA)l  =  (20  -  FOM)x  P{pfpf\p . . .) 

+  (FOTH)xP(Pfpf\pf..) 
+  (32  -  FOV)xP{pfpf\p.f.) 
+  (FOA)xP(pfpf\p..f) 

where  "  .  "  represents  either  a  pass  or  fail.    For  example, 

P(PPPP\P---)IS  tne  probability  that  device  passes  all  environments  given  that  it 
passed  manufacturer  test. 

P  {p  p  p  p\  p  p  •  .)  is  the  probability  that  a  device  passes  all  environments  given  that 
it  passed  manufacturer  test  and  temperature-humidity  environment. 

The  conditional  probabilities  are  computed  from  the  estimated  cell  probabilities 
from  the  previous  M  step.   As  an  example,  we  have  following  failure  vector, 

(1      1     2     0) 
Then  E  [  M      |  (  1,  1,  2,  0  )  ]  is  calculated  as  follows: 


'(«"!'->-   S   =0^67.93 


'(""I"")'     S    =0.9542328 


'(""I"-)-    gjggg    -0.9661516 


/>""'^l^")=    Si    °  0.937.388, 


where  the  numbers  come  from  (j.22)  the  estimated  cell  expected  values  after  from  the 
previous  M  step.   Thus, 
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ElMpppp\(l,  1,  2,  0)]  =  83.4654163 

The  remaining  16  conditional  expectations  are  calculated  similarly.  Then  we  com- 
pare these  expectations  with  previous  expectations.  If  the  difference  between  compared 
expectations  are  small  enough  then  EM  Algorithm  is  assumed  to  have  converged. 

When  EM  Algorithm  converges,  we  can  estimate  reliability  of  device  (  R  )  after 
exposure  all  four  environments  as 

A 

a  £1    Mpppp   \ 

R  =  920  •  (3-25) 

For  example,  for  the  failure  vector  above,  after  27  iterations,  we  obtain 

EL^pppp1  =  84.412714, 
and 


A  84  412714 

R  =         92  Q         =  °-9175259- 


D.     INITIAL  GUESS  PROBLEM 

The  EM  Algorithm  uses  an  initial  guess  vector  as  in  equation  (3.9).    We  use  the 
independence   assumption  as  described  in  equations  (3.1),  (3.2)  and  (3.3).   The  following 
procedures  are  used  during  the  application  of  EM  Algorithm  to  the  log  linear  model  with 
the  independence  assumption. 

We  have  still  an  initial  guesses  problem  for  parameters  in  the  log  linear  model  with 
the  independence  assumption.  For  this  reason,  we  do  some  calculations  to  find  initial 
guesses  for  the  parameters.  Moreover,  one  can  choose  random  initial  guesses  for  the 
parameters  too. 
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Assume  that  we  have  the  following  failure  vector, 

(10     10). 

In  this  procedure  we  divide  the  number  of  devices  which  pass  the  manufacturer  test 
by  the  number  of  cells  which  include  items  that  have  passed  the  manufacturer  test.  In 
this  example,  we  can  calculate  Mpp!f  as  follows: 

19  0 

•  [        '      ]  =  2.375  items  for  manufacturer  test 

O.O 

•  [        '      3  =  2.500  items  for  temperature  and  humidity  test 

O.O 

•  [_  -^—  ]  =  0.125  item  for  vibration  test 

O.O 

•  [  -tttt  ]  =  0.000  item  for  altitude  test 
finally  the  expected  number  of  device  for  this  cell  is, 


Mppff  =  5.000 


After  computing  M,jhl  for  each  cell  using  the  same  procedure  above,  we  can  group  these 
cells  to  estimate  Mpp00 ,  MfJM  ,  M/p00 ,  Mm  as  follows: 

•  YM„u  =  Mpp00 

•  YMfiki  =  Mm 

where  k,   /   =    p,  f. 
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In  this  example, 


Table  6.     GRUOPED  DATA  FOR  INITAL  GUESS  CALCULATION 


TERMS 

ESTIMATED  VALUES 

K^ 

32.500 

Kfio 

22.500 

Mfroo 

23.500 

Mff» 

13.500 

By  combining  the  parameters   n   and  X\ ,  we  have  four  equations  and  four  unknown 
parameters.   We  can  easily  solve  for  the  parameters  from  the  initial  guesses. 

Table  7.     ESTIMATED  PARAMETERS  FOR  INITIAL  GUESS 


TERMS 

ESTIMATED  VALUES 

A*  .+  K 

-1.97 

K 

1.08 

K 

1.08 

K 

0.56 

Finally,  we  can  split  the  sum  of  the  n  and  X\  into  two  parts.  One  possibility  is  to 
divide  the  sum  by  two  assigning  half  to  //  and  half  to  X\.  Results  using  a  random  initial 
guess  and  the  above  procedure  for  an  initial  guess  are  very  close.  Calculation  of  initial 
guess  parameters  are  done  by  the  FORTRAN  program  INITIAL  in  Appendix  D  and 
the  FORTRAN  program  PARAM  in  Appendix  E. 
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E.     RESULTS 

Reliabilities  of  the  pyrotechnic  device  are  calculated  by  a  FORTRAN  program 
LLMDEP  in  Appendix  F.  They  are  given  in  Table  8.  Results  from  the  log  linear  model 
design  with  dependence  assumption  are  similar  to  the  first  model.  We  can  easily  see 
from  the  previous  table  that  the  ordering  of  the  estimates  of  R  is  counter-intuitive.    We 

A  A 

expect  (  1  1  2  1  )  to  yield  the  smallest  R  rather  than  an  R  =  0.934  which  is  larger  than 
R  for  about  half  of  the  24  realizations  of  the  failure  vector. 

It  is  clear  from  these  results  that  this  log  linear  model  is  inappropriate  for  modeling 
the  outcomes  of  the  sampling  inspection  plan.  More  realistic  models  would  include 
three-way  and  four- way  interaction  terms.  However,  due  to  the  extreme  censoring  in  the 
data,  we  can  not  estimate  R  for  these  models.  Moreover,  reliabilities  of  several  cases  for 
even  this  model  were  not  calculated  because  of  computational  limitations.  Taking  the 
account  of  dependence  with  log  linear  models  is  clearly  not  a  reasonable  approach  for 
estimating  reliabilities  from  this  data. 
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Table  8.     RELIABILITIES  WITH  LOGLINEAR  MODEL 


FAILURE  VECTOR 

A 

(000  0) 

1.0000000 

(10  00) 

0.9887247 

(10  10) 

0.9693505 

(1001) 

0.9518477 

(1100) 

0.9503552 

(1110) 

0.9457148 

(10  11) 

0.9456805 

(1111) 

0.9399773 

(1020) 

0.9374950 

(112  1) 

0.9336501 

(0011) 

0.9255090 

(0  110) 

0.9255087 

(110  1) 

0.9246104 

(10  2  1) 

0.9175486 

(1120) 

0.9175295 

(0101) 

0.9089977 

(002  1) 

0.8978485 

(0  120) 

0.8978109 

(0111) 

0.8825021 

(0121) 

0.8584540 
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IV.     WORST  CASE  SCENARIO 

A.     ASSOCIATION  ANALYSIS 

In  the  inspection  sampling  plan,  failures  are  not  assigned  to  an  individual  failure 
mechanism  except  for  the  manufacturer  test.  For  example,  from  the  following  failure 
vector. 


(      1 


1 


1      ) 


We  know  that  there  is  at  least  one  failure  due  to  the  manufacturing  test.  There  is  one 
failure  from  the  joint  temperature  and  humidity  and  manufacturer  test,  but  we  do  not 
know  which  failure  mechanism  generated  this  failure.  Items  can  fail  due  to  a  manufac- 
turing related  failure  mechanism  or  other  failure  mechanism  or  both.  The  worst  case  is 
to  assume  that  the  cause  of  failure  is  due  to  all  of  the  failure  mechanisms  that  it  was 
exposed  to. 

Let  MT  represents  manufacturer  test,  TAHT  represents  temperature  and  humidity 
test,  VT  represents  vibration  test,  AT  represents  altitude  test. 

Because  all  92  items  are  exposed  to  manufacturer  test,  the  worst  case  of  the  failure 
analysis  of  the  failure  vector  above  is  shown  by  the  following  table. 

Table  9.     FAILURE  ANALYSIS 


MT 

TAHT 

VT 

AT 

#  OF  TESTED  ITEMS 

92 

20 

32 

20 

#  OF  SUCCESFUL 
ITEMS 

87 

19 

30 

19 

This  "  worst  case  "  scenario  should  give  us  lower  bounds  for  the  estimates  R  of  R. 
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If  we  calculate  the  reliability  of  series  system  assuming  failure  mechanisms  in  test 
are  independent  when  in  fact  they  are  associated  but  not  independent,  then  we  underes- 
timate system  reliability.  It  is  reasonable  to  expect  that  there  is  positive  dependence 
between  four  tests.  If  one  item  fails  in  manufacturer  test,  it  is  more  likely  to  fail  in  other 
tests.  One  way  of  modelling  positive  dependance  between  tests  is  to  assume  that  tests 
are  positively  quadrant  dependent. 

Given  random  variables  Tl  ,  . . . ,  Tn.  They  are  said  to  be  positively  quadrant  de- 
pendent (PQD)  [Ref.  8:  p.    33],  if 

n 

P(T,<  flf  ...,  r„  <  /„)  >  ]]/>(7}  <  tt)  (4.1) 

for  all    ( /,,  t2,  .  .  .  ,  /„)  e  R".    An  equivalent  formulation  of  positive  quadrant  depend- 
ence is  7i,  T2,  . . . ,  T"  are   PQD   iff 

n 

P(  Tx  >  r, ,  . . . ,  T„  >  tn)  2>  ]  ]/>(  Tt  >  rf)  (4.2) 

for  all    ( r„  t2,  . . . ,  /„)  e  Rn.    The  proof  that  (4.1)  and  (4.2)  are  equivalent  as  given  in 
[Ref.  8:  pp.    32-33]. 

We  may  take  account  positive  dependence  between  tests  with  the  PQD  assumption. 
Let 


{1  ,       if  an  item  passes  environment  i 
(4-3) 
0  ,       if  an  item  fails  environment  i 


for  i  =  1,2,3,4.  Assume  that  T, ,  T2 ,  T3 ,  T4  are  positively  quadrant  independent,  i.e. 
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P(  Tx  <  /,  ,  T2  <  t2 ,  T3  <  t3 ,  T4  <  t,)  >  ]  ]/>(  Tt  <  tt)  (4.4) 


j=i 


or  equivalently  that 


P(  T,  >  tx  ,  T2  >  t2  T3  >  t3 ,  T4  >  r4)  >  ]  "]/»(  7}  >  f,).  (4.5) 

Let  R  be  the  probability  that  an  item  activates  after  exposure  to  all  4  environments. 
Then 

R  =  />( tx  -  i ,  r2  =  i ,  r3  -  i ,  r4  =  1 ).  (4.6) 

Using  equations  (4.4)  and  (4.5) 

/?  >  p(  r,  =  1 )  p{  t2  =  i  )  />(  r3  -  i  )  P(  r4  =  1 ) .  (4.7) 

Using  the  notation  from  the  previous  section 

R  >  0i  Qi  Qi  Qa  ,  (4.8) 

where  QL  is  the  probability  that  an  item  passes  environment  i.  With  the  censored  data, 
we  can't  estimate  Q{ ,  Q2,  Q3  and  Q*  without  building  a  more  structured  model  for 
7", ,  T2 ,  T3  and  T4.  One  alternative  is  get  a  lower  bound  for  estimates  of 
Q\ »  Qi »  (23 .  Qa  •  This  is  worst  case  scenario  . 

B.     CALCULATIONS  WITH  EXAMPLE 

Let 

■  R  be  the  MLE  for  the  reliability  of  the  device. 

■  Qi  be  the  MLE  for  the  probability  that  an  item  passes  environment  i. 
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■  R  be  a  lower  bound  of  R  for  reliability  of  the  device. 

■  Q,  be  a  lower  bound  for  the  probability  that  an  item  passes  environment  i. 
Then  according  to  equation  (4.8) 


4 

*  >      \Qi 


i=\ 


Thus 


4 


A 
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(4.9) 


Further  we  will  construct  estimates  Q,  of  Q,  such  that 

&  <  Qt  (4.10) 

and  then  define 


R=         Qt.  (4.11) 


R  <  R.  (4.12) 


From  Chapter  II,  the  likelihood  of  observing  Xx  =  xu  X2  =  jc2,  X3  =  x3,  X4  =  x4  is 


L(x],  x2,  x2,  x4  |  R},  R2,  /?3,  /?4 )  — 


no** a -^ 


p{tx  =  i)Xi?(r,  =  o)' 


i-jii 


x  ("2JP(T]  =  l,T2=\)X2 

x  (1  -P{TX  =  \,T2=  l))"2'*2  (4-13) 

x  r3jp(r1  =  i,r3  =  i)X3 

x  (i  -/>(r1  =  i,r3  =  i))''3~JC3 
x  r4jp(r,  =  1,74=1)^ 
x  (i  -/>(r1  =  i,74=  l))**-* 

where  x,  is  the  number  of  devices  out  of  nl  that  activate  after  test  i.  If  we  know  why  the 
device  failed  for  the  tests  i  =  2,  3,  4  which  include  manufacturer  test  along  with  exposure 
to  environment  i  then  our  likelihood  could  be  written  as 
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L(X],  x2,  x3,  x<\Ru  R2,  R3,  P4)  =  )P(T]  =  l)x^P(T]=0) 


nx-xx 


\x2>   x2\'   •r22'    •r23' 


x  p (r,  =  i,  t2  =  \)*2 p (r,  =  i,  r2  =  0)*21 
x  />  (r,  =  o,  t2  =  if22  p  (r,  =  o,  r2  =  of23 

X  (  J  (4-14) 

\x3>   *31»   *32»   -*33/ 

x  p (Tj  =  i,  r3  =  if3  p (r,  =  i,  r3  =  of31 
x  p  (7j  =  o,  t3  =  i)X32  p  (r,  =  o,  r3  =  o)*33 

-(  "  ) 

\X4,    JT4j,   ^42,    -^43/ 

x  p (r,  =  i,  r4  =  i)*  p (r,  =  i,  r4  =  o)^1 
x  p  (r,  =  o,  r4  =  l)*42  p  (^  =  o,  r4  =  o)^3 

where  or,,  is  the  number  given  test  i  which  failed  due  to  environment  i  but  passed  manu- 
facturing, xa  is  the  number  given  test  i  that  passed  environment  i  but  failed  manufac- 
turing and  xa  is  the  number  given  test  i  that  failed  both  manufacturing  and  environment 
i  for  i  =  2,  3,  4.   Note  that  «,  —  x,  =  xn  +  xa  +  xa  for  i  =  2,  3,  4. 
From  this  likelihood,  the  MLE's  of  Q{  =  P{T{=  1)  are  given  by 


X,     +    (X2+X7A    +    (x-y  +  X-,,)    +    (xa  +  Xai) 

Q,  =  — — ^r — -r -r — — ,  (4.15) 
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Qi  =  n,         >  (4-16) 


~  X-i     ~T     X-i-) 


~  Xa      1      Xa-) 


However  since  we  do  not  know  the  je,/s,  we  see  that 


x,   +  x2  +  x2  +  x4 

ft  > -T-2 ~  (4-19) 


and 


ft  *    7T,  (4-20) 


for  i  =  2,  3,  4.    Let 


a  =  *■ +  *+"  +  *  (4.2.) 


and 


a-  =  it  (4-22) 

A  — 

for  i  =  2,  3,  4.   Then  <2/s  are  lower  bounds  for  the  true  MLE's  Q„  and 

*  =  ftftftft  (4-23) 

It  is  clear  that  R  is  a  lower  bound  for  the  MLE  R. 
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In  Worst  Case  Scenario,  we  assume  that  the  item,  which  fails  in  test  i,  fails  because 
of  both  manufacturing  related  failure  mechanism  and  other  failure  mechanism  which  is 
induced  by  the  i  th  test  environment.  With  this  assumption,  the  reliability  of  of  device 
is  given  as  follows 


A        A       A        /\ 


R  =  Q\  Qi  Q3  &  (4.24) 


where 


a  X\    +    X-)    +    X->    +    Xa 

Qx  =  -J— 


92 


(4-25) 


A  X-y 

A  JC, 

*"  If 

^4  "  20 


Here  is  an  example,  in  this  example  our  failure  vector  is, 

(0        1        2        1     ) 


&  =     20+19  +  30+19     =  0.9456521 

<?2  =   4tt   =  0.9500000 
^2         20 

4=4?-=  0.9375000 
J         32 

4  =   -^-   =  0.9500000 


And  finally  reliability  of  the  device  can  be  estimated  by  using  equation  (4.24) 

R  =  ( 0.9456521 )  ( 0.9500000 )  (  0.9375000 )  (  0.9500000 )  =  0.8093069      (4.27) 
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(4.26) 


C.     RESULTS 

We  compute  MLE's  with  a  FORTRAN  program  MLEB  in  Appendix  G.  These 
are  given  Tables  10  and  11. 

To  get  an  approximate  lower  confidence  bounds  for  R  using  the  worst  case  data, 
we  bootstrap  using  the  procedure  described  in  the  previous  chapter.  The  FORTRAN 
program  RANVEC  in  Appendix  B  is  used  to  generate  the  5000  bootstrap  samples  for 
each  case.  We  then  estimate  R  's  for  5000  failure  vectors  by  the  means  of  MLEB  for 
given  case.  The  next  step  is  to  compute  the  order  statistics  of  R  's  from  Rx  to  Rsm.  This 
is  done  by  FORTRAN  program  SORT  in  Appendix  C.  Finally,  we  obtain  the  95  % 
lower  confidence  bound  from  this  routine.  Reliabilities  and  95  %  lower  confidence 
bounds  are  tabulated  in  the  following  pages.  Results  of  reliabilities  and  95  %  lower 
confidence  bounds  are  given  in  descending  order. 

These  estimates  and  LCB's  for  R  are  decreasing  with  the  total  number  of  the  failures 
out  of  92  items  tested.  This,  at  least,  is  consistent  with  how  we  believe  pyrotechnic  de- 
vices behave.  It  should  be  noted  that  these  estimates  are  in  fact  conservative  lower 
bounds  for  the  time  MLE's  under  the  very  weak  assumption  of  PQD.  How  conservative 
these  estimates  are  cannot  be  determined  without  more  extensive  data  that  allows  us  to 
estimate  the  degree  of  dependence  between  tests. 
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Table  10.  TEST  PROBABILITIES  (  B  ) 

CASE 

K, 

k 

R, 

K 

(0  000) 

1.0000000 

1.0000000 

1.0000000 

1.0000000 

(1000) 

0.9891340 

1.0000000 

1.0000000 

1.0000000 

(0010) 

0.9891340 

1.0000000 

0.9687500 

0.1.00000 

(10  10) 

0.9782608 

1.0000000 

0.9687500 

1.0000000 

(0001) 

0.9891304 

1.0000000 

1.0000000 

0.9500000 

(0  100) 

0.9891340 

0.9500000 

1.0000000 

1.0000000 

(1001) 

0.9782608 

1.0000000 

1.0000000 

0.9500000 

(1100) 

0.9782608 

0.9500000 

1.0000000 

1.0000000 

(0  02  0) 

0.9782608 

1.0000000 

0.9375000 

1.0000000 

(1020) 

0.9673913 

1.0000000 

0.9375000 

1.0000000 

(00  11) 

0.9782608 

1.0000000 

0.9687500 

0.9500000 

(01 10) 

0.9782608 

0.9500000 

0.9687500 

1.0000000 

(10  11) 

0.9673913 

1.0000000 

0.9687500 

0.9500000 

(1110) 

0.9673913 

0.9500000 

0.9687500 

1.0000000 

(0101) 

0.9782608 

0.9500000 

1.0000000 

0.9500000 

(110  1) 

0.9673913 

0.9500000 

1.0000000 

0.9500000 

(002  1) 

0.9673913 

1.0000000 

0.9375000 

0.9500000 

(0  120) 

0.9673913 

0.9500000 

0.9375000 

1.0000000 

(1021) 

0.9565217 

1.0000000 

0.9375000 

0.9500000 

(1120) 

0.9565217 

0.9500000 

0.9375000 

1.0000000 

(0111) 

0.9673913 

0.9500000 

0.9687500 

0.9500000 

(1111) 

0.9565217 

0.9500000 

0.9687500 

0.9500000 

(0121) 

0.9565217 

0.9500000 

0.9375000 

0.9500000 

(112  1) 

0.9456521 

0.9500000 

0.9375000 

0.9500000 
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Table  11.     RELAIBILITIES     AND 
BOUNDS ( B ) 


95     %     LOWER     CONFIDENCE 


FAILURE 
VECTOR 

R 

FAILURE 
VECTOR 

95  %  LCB 

(0  000) 

1.0000000 

(00  00) 

1.0000000 

(1000) 

0.9891304 

(1000) 

0.9891304 

(00  10) 

0.9582200 

(00  10) 

0.8766983 

(10  10) 

0.9476901 

(10  10) 

0.8766983 

(0001) 

0.9396738 

(0001) 

0.8222825 

(0100) 

0.9396738 

(000  1) 

0.8222825 

(1001) 

0.9293478 

(0  100) 

0.8222825 

(1100) 

0.9293478 

(1  100) 

0.8222825 

(0020) 

0.9171195 

(0020) 

0.8179346 

(102  0) 

0.9069293 

(1020) 

0.8084239 

(0011) 

0.9003056 

(0011) 

0.7712974 

(01  10) 

0.9003056 

(0  110) 

0.7712974 

(10  11) 

0.8903022 

(10  11) 

0.7712974 

(1110) 

0.8903022 

(0101) 

0.7712974 

(0101) 

0.8828803 

(0101) 

0.7548369 

(110  1) 

0.8730705 

(110  1) 

0.7548369 

(0  120) 

0.8615828 

(0  120) 

0.7218070 

(002  1) 

0.8615828 

(002  1) 

0.7200747 

(1021) 

0.8519021 

(1021) 

0.7200747 

(1120) 

0.8519021 

(1120) 

0.7200747 

(0111) 

0.8457870 

(0111) 

0.6929346 

(1111) 

0.8362838 

(1111) 

0.6927614 

(0121) 

0.8093069 

(0121) 

0.6508338 

(112  1) 

0.8001103 

(112  1) 

0.6505433 
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V.     BONUS  SYSTEM  APPROACHES 

A.     BACKGROUND 

Quality  is  described  as  "  especially  high  degree  of  goodness  or  worth  "  [Ref.  9 
:  p.  685].  In  industry,  a  quality  product  is  one  that  fulfills  customer  expectations. 
There  are  two  general  aspects  of  quality: 

■  Quality  of  Design 

■  Quality  of  Conformance. 

All  goods  and  services  are  produced  in  various  grades  or  levels  of  quality.  These 
variations  in  grades  or  levels  of  quality  are  intentional;  therefore  the  appropriate  tech- 
nical term  is  quality  of  design.  The  quality  of  conformance  is  how  well  the  product 
conforms  to  the  specifications  and  tolerances  required  by  the  design.  Quality  of  con- 
formance is  influenced  by  a  number  of  the  following  factors: 

■  the  choice  of  manufacturing  process 

■  the  training  and  supervision  of  the  workforce 

■  type  of  quality  assurance  system  (  process  controls,  tests,  inspections,  etc.  )  used 

■  the  extent  to  which  these  quality-assurance  procedures  are  followed 

■  the  motivation  of  workforce  to  achieve  quality. 

Quality  Control  is  the  engineering  and  management  activity  by  which  we  measure 
the  quality  characteristics  of  a  product,  comparing  them  with  specifications  or  require- 
ments and  taking  appropriate  remedial  action  whenever  there  is  a  difference  between  the 
actual  performance  and  the  standard  [Ref.  10  :  pp.  1-3]. 

As  contracts  are  now  written,  contractors  need  only  to  satisfy  the  requirments  of  the 
sampling  inspection  plan  for  lot  acceptance.  Contractors  have  no  incentive  to  improve 
the  quality  of  the  items  they  provide,  although  they  are  in  a  position  to  do  so.    As 
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mentioned  before,  the  quality  of  conformance  is  influenced  by  the  motivation.  Therefore 
we  can  motivate  manufacturers  by  giving  a  bonus  for  improved  quality.  To  improve 
quality,  The  Naval  Weapons  Support  Center  has  decided  to  implement  a  bonus  system 
for  pyrotechnic  devices. 

B.     BONUS  PLANS 

In  this  chapter,  we  design  a  bonus  system  to  improve  the  quality  of  pyrotechnic 
devices.  A  good  bonus  system  encourages  reliability  growth,  because  firms  try  to  reach 
a  high  quality  to  get  a  bonus.  An  effective  bonus  system  must  detect  small  differences 
among  the  offered  lots.  In  the  previous  chapters,  we  estimated  the  reliability  of  a 
pyrotechnic  device  in  several  ways.  We  assumed  independence  in  the  first  model,  and 
we  assumed  dependence  between  tests  in  the  last  two  models.  The  estimated  reliabilities 
in  the  first  two  models  are  close  to  each  other,  but  they  exhibit  different  structures  in 
order.  One  can  easily  see  that  the  failure  vector  of  the  worst  case  (1  12  1)  has  bigger 
reliability  value  than  ten  of  the  possible  cases,  from  Table  5  and  Table  15.  In  addition, 
the  order  of  95  %  LCB's  of  cases  does  not  match  to  the  order  of  MLE's  in  independent 
models.  But  in  the  worst  case  scenario,  which  assumes  dependence  between  tests,  we 
get  the  same  order  for  both  95  %  LCB's  and  MLE's.  Thus,  we  will  use  the  worst  case 
scenario  model  in  further  calculations  for  the  bonus  system.  With  the  assumptions 
above,  we  can  apply  three  sampling  plans  for  giving  bonuses  to  the  manufacturers. 
These  are: 

■  Single  Sampling  Bonus  System 

■  Double  Sampling  Bonus  System 

■  Multi-Sampling  Bonus  System. 

Manufacturers  have  to  meet  pyrotechnic  device  acceptance  criteria  first,  before 
having  a  chance  to  get  a  bonus. 
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1.     Single  Sampling  Bonus  System 

In  this  bonus  system,  we  will  decide  two  things  at  the  end  of  inspection.  First 
we  will  decide  whether  the  ofTered  lot  is  acceptable  or  not,  and  then  for  lots  which  are 
accepted  we  will  decide  whether  to  give  a  bonus  to  the  manufacturer. 

The  Single  sampling  bonus  system  is  designated  by  three  numbers.   These  are, 

«,  ,  LCBFB  ,  LCB 

where 

■  «!  is  the  sample  size  for  inspection  (  here  «,  =  92  ). 

■  LCBFB  is  the  cut  off  value  applied  to  the  lower  confidence  bound  for 
awarding  a  bonus. 

■  LCB  means  estimated  lower  confidence  bound  after  inspection.    (If  LCB  is 
greater  than  or  equal  to  LCBFB,  then  a  bonus  is  awarded.) 

The  following  algorithm  shows  us  how  single  sampling  bonus  system  works  for 
pyrotechnic  devices. 

STEP  #  1   :   DETERMINE  CUT  OFF  VALUE  FOR  AWARDING  THE  BONUS. 
STEP  #  2   :   TAKE  A  SAMPLE  SIZE  OF  92  FROM  OFFERED  LOT. 
STEP  #  3   :   APPLY  MANUFACTURER  AND  THREE  ENVIRONMENT  TESTS. 
STEP  #  4   :   COMPARE  RESULTS  OF  TESTS  WITH  ACCEPTANCE  CRITERIA. 
STEP  #  5   :   ESTIMATE  ITS  LOWER  CONFIDENCE  BOUND  (  LCB  ) , 

IF  THEY  MET  ACCEPTANCE  CRITERIA. 
STEP  #  6   :   COMPARE  LCB  OF  LOT  WITH  LCBFB 
STEP  #7   :   IF  THE  LCB  IS  GREATER  THAN  OR  EQUAL  TO  LCBFB, 

GIVE  BONUS  TO  THE  FIRM. 
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According  to  the  algorithm  described  above,  the  following  are  possible  events 
for  the  firms. 

■  Firm  may  not  satisfy  our  acceptance  criteria.    This  means  that  firm  gets  a 
failure  vector  worse  than   (112  1). 

■  Firm  satisfies  acceptance  criteria,  but  its  LCB  may  be  less  than  LCBFB.   This 
means  that  firm  does  not  get  bonus,  but  the  lot  is  accepted. 

■  Firm  satisfies  acceptance  criteria,  and  its  LCB  may  be  greater  than  or  equal 
to  LCBFB.   This  means  that  firm  gets  the  bonus. 

2.     Double  Sampling  Bonus  System 

A  double  sampling  plan  has  an  advantage  over  a  single  sampling  plan.  Because 
a  double  sampling  plan  involves  a  larger  sample  size,  it  reduces  the  chance  that  a  man- 
ufacturer who  deserves  a  bonus  will  not  get  one.  Double  Sampling  Bonus  System  per- 
mits the  taking  of  two  samples  on  which  to  make  a  decision  [Ref.     11  :  pp.   184-185]. 

In  this  system,  we  have  two  inspection  stages.  If  the  firm  does  not  get  a  bonus 
after  the  first  inspection,  then  the  firm  is  given  a  second  chance  with  a  second  inspection. 
A  double  sampling  bonus  system  is  designated  by  five  numbers. 

«,  ,  n2 ,  LCBFB ,  LCBl  ,  LCB2 

where 

■  «,  is  the  sample  size  for  first  inspection,  (  here  nx  =  92  ) 

■  «2  is  the  sample  size  for  second  inspection,  (  here  n2  =  92  ) 

■  LCBFB  is  the  cut  off  value  applied  to  lower  confidence  bounds  for  awarding 
bonuses, 

■  LCBl  is  the  estimated  lower  confidence  bound  after  first  inspection, 

■  LCB2  is  the  updated  estimate  of  lower  confidence  bound  after  second  in- 
spection. 
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In  the  double  sampling  bonus  system,  if  the  firm  does  not  get  a  bonus  after  the 
first  inspection  but  does  meet  the  acceptance  criteria  for  the  first  sample,  then  a  second 
sample  is  taken.  We  already  have  LCB's  of  possible  cases  after  the  first  inspection  in 
Table  1 1.  After  the  second  inspection,  we  calculate  the  LCB  using  an  aggregated  failure 
vector  which  includes  failures  from  both  samples.  When  we  compute  LCB's  using  ag- 
gregated failure  vector  after  the  second  inspection,  we  have  a  total  sample  size  of  184; 
there  are  135  different  failure  vectors  for  which  the  lot  meets  the  acceptance  criteria  for 
both  samples.  After  tabulating  these  possible  135  cases  and  the  estimates  of  R,  we  used 
the  bootstrap  procedure  to  find  LCB's  for  each  case.  We  created  5000  random  failure 
vectors  for  each  of  them  by  using  case  success  probabilities  in  tests.  After  this,  we  esti- 
mated MLE's  of  5000  failure  vectors  for  each  possible  case.  Finally,  we  estimated  95 
%  LCB  of  each  case.  These  calculations  were  done  by  using  programs  RANVEC, 
MLEB,  SORT  in  appendix  B,  G  and  C  respectively.  The  results  are  tabulated  at  Ap- 
pendix I.  The  following  algorithm  shows  us  how  double  sampling  bonus  system  works 
for  pyrotechnic  device. 

STEP  #      1    :      DETERMINE  LCBFB. 

STEP  #     2    :      TAKE  A  SAMPLE  SIZE  OF  92  FROM  OFFERED  LOT. 

STEP  #     3    :      APPLY  MANUFACTURER  AND  THREE  ENVIRONMENT  TESTS. 

STEP  #     4    :      COMPARE  RESULTS  OF  TESTS  WITH  ACCEPTANCE  CRITERIA. 

STEP  #     5    :      ESTIMATE   ITS  LOWER  CONFIDENCE   BOUND   (   LCB1    ), 

IF  THEY  MET  ACCEPTANCE  CRITERIA. 
STEP  //     6    :      COMPARE  LCB1  OF  LOT  WITH  LCBFB 
STEP  #     7    :      IF  THE  LCB1   IS  GREATER  THAN  OR  EQUAL  TO  LCBFB, 

GIVE  BONUS  TO  THE  FIRM. 
STEP  #     8    :      IF  THE  LCB1   IS  LESS  THAN  LCBFB, 
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GIVE  A  SECOND  CHANCE  TO  THE  FIRM  FOR   BONUS. 
STEP  #     9    :      TAKE  A  NEW  SAMPLE  SIZE  OF  92  FROM  OFFERED  LOT. 
STEP  #    10    :      APPLY  MANUFACTURER  AND   THREE   ENVIRONMENT  TESTS 
STEP  #   11    :      AFTER  GETTING  THE  NEW  FAILURE  VECTOR, 

ADD  THIS  ONE  TO  THE  FIRST  FAILURE  VECTOR. 
STEP  #    12    :      ESTIMATE    ITS  LOWER  CONFIDENCE   BOUND   (    LCB2    ), 

WITH  AGGREGATED  FAILURE  VECTOR. 
STEP  #    13    :      COMPARE   LCB2   OF  LOT  WITH  LCBFB 
STEP  #    14    :      IF  THE   LCB2    IS  GREATER   THAN   OR  EQUAL  TO  LCBFB, 

GIVE   BONUS   TO  THE  FIRM. 


According  to  the  algorithm  described  above,  the  following  are  possible  events 
for  the  firms. 

■  The  firm  may  not  satisfy  our  acceptance  criteria  after  first  inspection.  This 
means  that  firm  gets  a  failure  vector  worse  than  (112  1).  The  lot  is  not  accepted  and 
second  sample  is  not  taken. 

■  The  firm  satisfies  acceptance  criteria,  and  its  LCB1  may  be  greater  than  or 
equal  to  LCBFB.  This  means  that  it  gets  bonus  after  first  inspection,  and  a  second 
sample  is  not  needed. 

■  The  firm  satisfies  acceptance  criteria,  but  its  LCB1  may  be  less  than  LCBFB. 
This  means  that  it  does  not  get  the  bonus  after  first  inspection,  but  still  has  a  chance  to 
get  a  bonus  if  it  submits  a  second  sample. 

Of  those  which  submit  a  second  sample 

■  LCB2  may  be  less  than  LCBFB.  It  means  that  firm  does  not  get  the  bonus, 
but  the  lot  is  accepted. 
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■  The  aggregated  failure  vector  meets  acceptance  criteria.    But  LCB2  may  be 
grater  than  or  equal  to  LCBFB.    It  means  that  firm  gets  bonus. 

3.     Multi-Sampling  Bonus  System 

In  this  system  approach,  we  may  take  a  sample  of  size  92  from  offered  lot  three 
times,  four  times  or  more.  But  if  we  use  the  three  sampling  bonus  system,  after  the  third 
inspection,  we  may  get  (4X4X7x4  =  448  )  448  different  failure  vectors.  Multi- 
Sampling  Bonus  System  is  going  to  be  more  computationally  intensive  with  respect  to 
double  sampling  bonus  system.  This  can  be  done  later  using  with  the  same  reasoning 
in  double  sampling  bonus  system. 

C.     EXAMPLES 

1.     Single  Sampling  Bonus  System 

Let  us  assume  that  we  decided  that  lower  confidence  bound  for  bonus  will  be 
0.8500.   We  are  going  to  use  this  LCBFB  in  these  examples. 

■  Case  #  1 

At  the  end  of  inspection,  firm  has  following  failure  vector 

(200      0    ). 

The  firm  did  not  meet  the  acceptance  criteria.    We  immediately  reject  the  lot  and  no 
bonus  is  given. 

■  Case  #  2 

At  the  end  of  inspection,  firm  has  following  failure  vector 

(10    0      1    ). 

The  firm  meets  the  acceptance  criteria.    We  estimate  its  95  %  lower  confidence  bound 
to  be 
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LCB  =  0.8222825  . 

Because  LCB  <  LCBFB  ,  we  do  not  give  a  bonus  to  the  firm. 

■  Case  #  3 

At  the  end  of  inspection,  firm  has  following  failure  vector 

(001      0    ). 

The  firm  meets  the  acceptance  criteria.    We  estimate  its  95  %  lower  confidence  bound 
to  be 

LCB  =  0.8766983    . 

Because  LCB  >  LCBFB  ,  we  give  a  bonus  to  the  firm. 
2.     Double  Sampling  Bonus  System 

■  Case  #  1 

At  the  end  of  first  inspection,  firm  has  following  failure  vector. 

(10    0      1    ). 

The  firm  meets  the  acceptance  criteria.    We  estimate  its  95  %  lower  confidence  bound 
to  be 

LCB\  =  0.8222825    . 

Because   LCBl  <  LCBFB  ,  we  give  a  second  chance  to  the  firm.    Here  is  the  result  of 
second  inspection 

(032      1    ). 
The  aggregated  failure  vector  will  be 
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(13     2      2    ). 

It  is  clear  that  LCB2  <  LCB1  for  this  data,  thus  we  do  not  give  a  bonus  to  the  firm. 

■  Case  #  2 

At  the  end  of  first  inspection,  firm  has  following  failure  vector 

(10    0      1    ). 

The  firm  meets  the  acceptance  criteria.    We  estimate  its  95  %  lower  confidence  bound 
to  be 

LCBX  =  0.8222825  . 

Here,   LCB  <  LCBFB  thus  we  will  give  a  second  chance  to  the  firm.    Here  is  the  result 
of  second  inspection 

(001      1    ). 

The  aggregated  failure  vector  will  be 

(10     1      2    ), 

and  we  estimate  its  95  %  lower  confidence  bound  to  be 

LCBX  =  0.8369564    . 

Because  LCB  <  LCBFB  ,  we  do  not  give  a  bonus  to  the  firm. 

■  Case  #  3 

At  the  end  of  first  inspection,  firm  has  following  failure  vector 

(10    0      1    ). 
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The  firm  meets  the  acceptance  criteria.    We  estimate  its  95  %  lower  confidence  bound 
to  be 

LCBl  =  0.8222825    . 

Because  LCBl  <  LCBFB  ,  we  will  give  a  second  chance  to  the  firm.    Here  is  the  result 
of  second  inspection 

(000      1    ). 

The  aggregated  failure  vector  will  be 

(10     0      2    ), 

and  we  estimate  its  95  %  lower  confidence  bound  to  be 

LCBl  =  0.8706521  . 
Because  LCBl  >  LCBFB  ,  we  will  give  a  bonus  to  the  firm. 
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VI.     SIMULATION  RESULTS  OF  BONUS  SYSTEMS 

A.     BACK  GROUND 

We  simulated  the  bonus  system  to  get  an  idea  of  how  well  the  bonus  system  pro- 
posed in  the  previous  chapter  works.  This  provides  the  user  with  a  means  of  setting  the 
cut  off  criteria  for  awarding  a  bonus.  To  generate  random  failure  vectors,  we  need  to 
know  the  probabilities  of  being  successful  in  each  test  for  the  firm.  First,  we  assume  that 
the  firm  has  an  equal  probability  of  being  successful  in  each  test,  with  the  following 
values  of  Px  =  P2  =  P3  =  P, 

•  0.9375 

•  0.9500 

•  0.9750 

•  0.9900 

•  0.9950. 

In  a  second  set  of  simulations,  we  assume  that  the  manufacturer  test  has  a  bigger 
probability  of  being  successful  than  the  other  environment  tests.  The  following  proba- 
bilities of  being  successful  in  manufacturer  test  were  used. 

•  0.9990 

•  0.9950 

•  0.9750 

•  0.9500. 
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In  these  cases,  we  used  equal  probabilities  of  being  successful  in  other  environment  tests 
with  values  between  0.9375  and  manufacturer  test  probability  for  that  given  case. 
Finally,  we  assume  that  the  manufacturer  test  should  have  smaller  probability  of  being 
successful.  When  we  study  these  cases  with  this  assumption,  we  used  probabilities  of 
being  successful  in  other  environment  tests  between  manufacturer  test  probability  and 
0.9990  for  that  given  case.  We  assume  that  minimum  value  of  the  manufacturer  test 
probability  will  be  0.9500  because  of  worst  case  in  last  two  assumptions. 

We  generated  2000  random  failure  vectors  for  each  possible  combination  of  proba- 
bilities by  using  program  RANVEC  in  Appendix  B.  They  were  used  for  1000  repli- 
cations of  each  bonus  system,  because  the  Double  Sampling  Bonus  System  (  DSBS  ) 
potentially  uses  two  failure  vectors  per  replication. 

After  getting  the  failure  vectors,  the  next  step  is  to  decide  lower  confidence  bound 
for  giving  bonus  (  LCBFB  ).  We  chose  0.800,  0.825,  0.850,  0.875,  0.900,  0.950,  0.999 
as  LCBFB  during  our  simulations.  We  simulated  bonus  systems  by  using  program  BO- 
NUS in  Appendix  G.  Program  BONUS  counts  how  many  times  firm  gets  bonus  during 
1000  replications.  And  then  it  calculates  the  bonus  percentage  dividing  counted  number 
by  1000.  For  each  scenario  this  bonus  percent  is  an  estimate  of  the  probability  of  getting 
a  bonus. 

B.     INITIAL  COMPARISON  OF  SYSTEMS 

In  this  part,  with  equal  probabilities  of  being  successful  in  tests,  we  tried  to  see  the 
difference  between  the  Single  and  the  Double  Sampling  Bonus  Systems.  For  this  reason, 
we  simulated  Single  Sampling  Bonus  System.  Results  were  tabulated  and  plotted  in  next 
pages. 
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Table   12.     SINGLE  SAMPLING  BONUS  SYSTEM  (EQUAL  PROBABILITIES) 


PROB.'S 

CHOSEN  LOWER  CONFIDENCE  BOUNDS  FOR  BONUS 

0.800 

0.825 

0.850 

0.875 

0.900 

0.950 

0.999 

0.9200 

0.025 

0.007 

0.007 

0.007 

0.002 

0.002 

0.002 

0.9375 

0.060 

0.023 

0.023 

0.023 

0.010 

0.010 

0.007 

0.9450 

0.092 

0.044 

0.044 

0.044 

0.018 

0.018 

0.010 

0.9500 

0.117 

0.053 

0.053 

0.053 

0.022 

0.022 

0.011 

0.9600 

0.222 

0.112 

0.112 

0.112 

0.045 

0.045 

0.028 

0.9700 

0.380 

0.198 

0.198 

0.198 

0.107 

0.107 

0.067 

0.9750 

0.482 

0.281 

0.281 

0.281 

0.163 

0.163 

0.067 

0.9800 

0.609 

0.381 

0.381 

0.381 

0.237 

0.237 

0.106 

0.9850 

0.740 

0.511 

0.511 

0.511 

0.362 

0.362 

0.166 

0.9900 

0.868 

0.646 

0.646 

0.646 

0.502 

0.502 

0.274 

0.9950 

0.955 

0.805 

0.805 

0.805 

0.502 

0.502 

0.407 

It  is  obvious  that  there  is  no  difference  between  some  lower  confidence  bounds  for  bonus 
from  the  table  above.  The  firm  gets  the  same  bonus  percentage  when  we  use  0.825, 
0.850,  0.875  as  lower  confidence  bound  for  bonus.  The  same  thing  occurs  when  we  use 
0.900  and  0.950  as  LCBFB.  For  this  reason,  we  are  going  to  see  four  curves  in  Figure 
1.    Figure  2   shows  Double  Sampling  Bonus  System  with  different  LCBFB's. 
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SINGLE  SAMPLING  BONUS  SYSTEM 
DIFFERENT  LCB'S  FOR  BONUS  WfTH  EQUAL  TEST  PROBABIUT1ES 
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Figure  1.      Single  Sampling  Bonus  System  With  Different  LCBFB's 


DOUBLE  SAMPLING  BONUS  SYSTEM 

DIFFERENT  LCB'S  FOR  BONUS  wrm  EQUAL  TEST  PR0BABIUT1ES 
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Figure  2.      Double  Sampling  Bonus  System  With  Different  LCBFB's 
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As  you  see  in  figure  2,  we  have  six  different  curves  for  Double  Sampling  Bonus  System. 
Because  in  this  system,  we  can  see  from  the  table  that  LCBFB's  0.900  and  0.950  have 
approximately  the  same  bonus  percentage.    For  this  reason  we  did  not  plot  for  0.950. 

Table   13.     DOUBLE  SAMPLING  BONUS  SYSTEM  (EQUAL  PROBABILITIES) 


PROB/S 

CHOSEN  LOWER  CONFIDENCE  BOUNDS  FOR  BONUS 

0.800 

0.825 

0.850 

0.875 

0.900 

0.950 

0.999 

0.9200 

0.026 

0.008 

0.007 

0.007 

0.002 

0.002 

0.002 

0.9375 

0.064 

0.030 

0.025 

0.023 

0.010 

0.010 

0.007 

0.9450 

0.108 

0.059 

0.050 

0.045 

0.020 

0.018 

0.010 

0.9500 

0.142 

0.079 

0.063 

0.057 

0.027 

0.022 

0.011 

0.9600 

0.277 

0.172 

0.144 

0.122 

0.052 

0.045 

0.028 

0.9700 

0.498 

0.351 

0.293 

0.236 

0.138 

0.107 

0.067 

0.9750 

0.624 

0.487 

0.423 

0.353 

0.212 

0.163 

0.106 

0.9800 

0.752 

0.632 

0.565 

0.487 

0.319 

0.237 

0.166 

0.9850 

0.861 

0.799 

0.738 

0.643 

0.486 

0.362 

0.274 

0.9900 

0.948 

0.920 

0.887 

0.799 

0.671 

0.502 

0.407 

0.9950 

0.987 

0.982 

0.970 

0.941 

0.872 

0.695 

0.626 

Figure  2  shows  DSBS  to  be  more  sensitive  in  the  sense  that  a  higher  percentage  of 
the  firms  were  awarded  a  bonus.  For  this  reason,  we  decided  to  implement  DSBS  in  all 
simulations. 
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C.     SIMULATION  RESULTS  WITH  DIFFERENT  LCB'S  FOR  BONUS 

In  this  section,  results  for  each  chosen  lower  confidence  bound  for  bonus  will  be 
presented  as  follows.  We  used  values  in  Table  13  to  draw  plots  with  equal  probabilities 
in  each  test.    Polynomial  approximation  was  used  in  curve  fitting. 


DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  ■=  0.800 
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Figure  3.      Double  Sampling  Bonus  System  With  LCBFB  =  0.800 


Table   14.     DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.800 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.800 

PROB/S 

BONUS  % 

PROB/S 

BONUS  % 

PROB/S 

BONUS  % 

0.9200 

0.026 

0.9375 

0.064 

0.9450 

0.108 

0.9500 

0.142 

0.9600 

0.277 

0.9700 

0.498 

0.9750 

0.624 

0.9800 

0.752 

0.9850 

0.861 

0.9900 

0.948 

0.9950 

0.987 
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In  Figure  3,  when  the  probabilities  of  being  successful  in  each  test  increase,  then  the 
bonus  percent  increases.  For  example,  when  the  test  probabilities  is  equal  to  0.9200, 
then  the  bonus  percent  is  0.026.  If  the  firm  increases  its  probabilities  of  being  succesful 
in  each  test  to  0.9900,  then  the  bonus  percent  becomes  0.948.  The  result  of  the  double 
sampling  bonus  system  with  different  lower  confidence  bounds  are  tabulated  and  plotted 
in  Appendix  J. 

To  see  whether  the  probability  of  passing  the  manufacturer  test  effects  the  bonus 
percentage  differently  than  the  other  probabilities,  we  assumed  that  the  environment 
tests  would  have  equal  probabilities.  But  the  firm  will  have  a  different  probability  of 
being  successful  in  the  manufacturer  test.  In  the  following  table,  the  first  row  represents 
probability  of  being  successful  in  the  environment  tests  and  the  first  column  represents 
probability  of  being  successful  in  the  manufacturer  test.  The  intersection  of  rows  and 
columns  gives  us  the  bonus  percentage  of  a  firm  with  given  probabilities  of  being  suc- 
cessful in  the  tests.  These  procedures  were  done  for  each  LCBFB  value  separately.  The 
results  are  tabulated  and  plotted  in  Appendix  K.  Bonus  percentages  for  LCBFB  = 
0.8000  are  plotted  on  the  following  page.  It  is  clear  that  the  bonus  percentage  of  the 
firm  will  be  high  if  the  firm  has  big  probabilities  of  being  succesful  in  both  manufacturer 
test  and  other  joint  environment  tests. 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 


LCBFB  =  0.800 


DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  =  0.800 


Figure  4.      Double  Sampling  Bonus  System  With  LCBFB  =  0.800 


Table   15.     DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.800 


0.950 

0.975 

0.990 

0.995 

0.950 

0.142 

0.476 

0.683 

0.700 

0.975 

0.176 

0.624 

0.869 

0.891 

0.990 

0.202 

0.688 

0.948 

0.973 

0.995 

0.206 

0.698 

0.962 

0.987 
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D.     BONUS  PERCENTAGE   (  BPRCT  )   FORMULATION 

In  this  section  we  try  to  approximate  probability  of  getting  a  bonus  as  a  function 
of  lower  confidence  bound  for  bonus  (  LCBFB  ),  probability  of  passing  manufacturer 
test  (  P\  )  ,  and  probability  of  passing  environment  tests  (  P2  ).  To  do  this,  we  use  re- 
gression analysis  with  GRAFSTAT  using  the  bonus  percentages  which  are  the  simulated 
values  of  the  probability  of  getting  a  bonus.  If  a  reasonable  relationship  between  the 
bonus  percentage  and  LCBFB  and  the  probabilities  /*,  and  P2  can  be  found;  it  can  be 
used  to  set  LCBFB  without  resorting  to  simulation. 

We  want  to  formulate  bonus  percentages  as  a  function  of  the  following: 

■  Lower  confidence  bound  for  bonus  (  LCBFB  ) 

■  Probability  of  passing  manufacturer  test  Pi 

■  Probability  of  passing  environment  tests  P2. 

After  polynomial  approximation,  we  can  see  that  plotted  graphs  (Figures  3-10)  look 
like  logistic  growth  curves  [Ref.  12:  p.  383].  Then  we  can  formulate  bonus  probabilities 
(BP)  using  the  logistic  growth  function. 


BP  =   !— r  (7.4) 

1  +  eA 


where 


A  =  p0  +  /?,  Px  +  P2  P\  +  (}3  P2  +  /?4  P22  +  P5  />,  P2  +  P6  LCBFB        (7.5) 


We  can  make  a  transformation  described  as  below. 


1  +  eA  =   -^  (7.6) 


a         1  —  BP  ,n  -x 

e     =   —sp—  (7.7) 
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A  - iog  I.  ^sF-  ]• 


(7.8) 


Now  we  have  a  linear  equation  as  a  function  of  LCBFB,  Px  and  P2.    We  may  use 
simulation  results  and  approximate  parameters  doing  a  linear  regression. 

■  Let  Y  be  the  (  n  x  1)  column  vector  of  observations  on  dependent  variable. 

■  Let  X  be  the  (  n  x  p'  )  matrix  consisting  of  a  column  of  ones,  which  is  labeled  1, 
followed  by  p  column  vectors  of  observations  on  independent  variables. 

■  Let  P  be  the  (  p'  x  1  )  vector  of  parameters  to  be  estimated. 

■  Let  £  be  the  (  n  x  1  )  vector  of  random  errors.   Then 

Y  =  Xp  +  c  (7.9) 

We  can  obtain  observations  on  dependent  variable  as  below, 

1  -  BPRCT   1 


=    log  [ 


BPRCT 


(7.10) 


where  BPRCT  is  the  bonus  percentage  obtained  from  the  simulation.    As  an  example 
we  have  the  following  information  from  simulations. 
.  LCBFB  =  0.875 
.  />,  =  0.9500 
.  P\  =  0.9025 
.  P2  =  0.9375 
.  P\  =  0.8789 
.  PiP2  =  0.8910 
.  BPRCT  =  0.029 
.  Y  =  3.5110. 
We  can  write  following  equation. 
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3.5110  =  p0  +  /?,  0.95  +  fi2  0.9025  +  /?3  0.9375  +  /?4  0.8789  +  /?5  0.8910  +  P6  0.875  +  £ 

We  deliberately  chose  30  random  results  from  simulations.  We  formulate  them  as 
above.  We  can  estimate  unknown  parameters  on  these  equations  doing  linear  re- 
gression.   30  equations  can  be  written  with  matrix  notation  as  follows: 


Y  = 


Y(l) 
Y{2) 
F(3) 

?'(30) 


(7.11) 


And  X  matrix  is  designed  as  follows: 


X  = 


1  P,(l)  tf(l)  P2(l)  Pfo)  PlP1{\) 
,  Pi(2)  P?(2)  /»2(2)  P22(2)  /»,P2(2) 
1    />,(3)     pfo)     P2(3)     p2(3)     />,/>2(3) 


1  P,(30)  />J(30)  />2(30)  P22(30)  />,  P2(30) 


(7.12) 


And  parameter  vector  will  be  as  below. 


P  = 


Po 

Pi 

Pi 

P3 

P. 

Ps 

Pe 

(7.13) 
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We  estimate  parameters  with  linear  regression  using  GRAFSTAT  packages.  And 
then  we  formulate  the  bonus  percentage  using  the  estimates  of  parameters.  The  esti- 
mated bonus  probabilities  are  given  by 


BP  = 


1 


1   +  e 


xp 


(7.14) 


The  following  estimates  of  the  parameters  were  obtained  from  the  regression. 


Table  16.     COEFFICIENT  VECTORS  OF  REGRESSION  ANALYSIS 


PARAMETERS 

ESTIMATED  VALUE 

Po 

-1042.3 

fit 

695.11 

A 

-86.661 

A 

1560 

A 

-569.81 

A 

-570.62 

A 

12.266 

The  standard  error  was  0.45471  after  the  linear  regression.  We  can  use  these  esti- 
mates of  parameters  in  equation  (7.5).  And  we  can  estimate  the  bonus  probability  of  the 
firm  when  the  probability  of  passing  tests  and  lower  confidence  for  bonus  are  known. 
Some  examples  are  calculated  with  equation  (7.5).  Results  of  the  calculations  and 
comparison  with  simulation  are  summarized  in  following  table. 
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Table   17.     BONUS  PERCENTAGES  WITH  REGRESSION  ANALYSIS 

LCBFB 

P> 

P2 

SIMU- 
LATION 

EST.VALUE 

0.800 

0.9750 

0.9750 

0.624 

0.564 

0.800 

0.9990 

0.9375 

0.106 

0.062 

0.825 

0.9750 

0.9950 

0.882 

0.911 

0.825 

0.8500 

0.9750 

0.064 

0.077 

0.850 

0.9500 

0.9900 

0.612 

0.615 

0.850 

0.9900 

0.9750 

0.476 

0.529 

0.875 

0.9000 

0.9750 

0.114 

0.080 

0.875 

0.9700 

0.9700 

0.236 

0.211 

0.900 

0.9950 

0.9990 

0.985 

0.945 

0.900 

0.8750 

0.9900 

0.158 

0.096 

0.950 

0.9450 

0.9450 

0.018 

0.010 

0.950 

0.9750 

0.9990 

0.839 

0.790 

0.999 

0.9850 

0.9850 

0.274 

0.335 

0.999 

0.9900 

0.9250 

0.004 

0.002 
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VII.     CONCLUSIONS  AND  RECOMMENDATIONS 

We  estimated  the  reliability  of  pyrotechnic  device  from  the  sampling  plan  using  the 
following  models: 

■  Maximum  likelihood  assuming  independence  of  pyrotechnic  device  activation  in 
different  environments; 

■  Log  linear  model  incorparating  some  dependence  between  pyrotechnic  device  ac- 
tivation in  different  environments; 

■  Worst  case  scenario  which  gives  a  lower  bound  for  the  estimated  reliability  for  the 
general  model  where  no  assumptions  are  made  about  the  form  of  dependence  between 
pyrotechnic  device  activation  in  different  environments. 

Using  these  models  and  based  on  sampling  plan  data,  estimates  of  overall  reliability 
along  with  95  %  lower  confidence  bound  are  obtained.  We  computed  the  lower  confi- 
dence bound  for  each  possible  case  by  bootstrapping. 

Results  from  the  first  model  are  not  consistent  with  the  way  pyrotechnic  devices 
operate.  In  particular,  the  estimated  reliabilities  are  not  ordered  as  we  expected  them 
to  be.  Intiutively,  we  expect  samples  with  fewer  failures  to  give  smaller  R  values  and 
corresponding  LCB's  than  samples  with  more  failures.  This  is  not  the  case  for  the  first 
model.  For  example,  the  failure  vector  that  has  the  maximum  number  of  failures  for 
each  test  (  a  total  of  5  failures  )  is  in  the  middle  of  the  order  with  respect  to  R  and  the 
95  %  lower  confidence  bound.  The  failure  vector  (0101)  with  a  total  of  2  failures 
has  a  much  lower  R  and  LCB.  The  discrepency  between  the  results  of  this  model  and 
what  we  expected  to  see  are  probably  due  to  the  fact  that  there  is  a  dependence  between 
the  events  which  a  device  activates  under  different  environments. 
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The  results  of  the  log  linear  model  design  with  two-way  dependence  assumptions  are 
similar  to  the  first  model.  Two-way  interaction  terms  were  used  in  this  model.  More- 
over, reliabilities  of  several  cases  were  not  calculated  because  of  computational  limita- 
tions. Due  to  the  extensive  censoring  in  the  sample  data,  it  does  not  appear  to  be 
possible  to  estimate  the  reliability  based  on  models  which  incorparate  dependence. 
Thus,  we  turn  to  finding  lower  bounds  for  the  estimated  reliability  based  on  models  with 
dependence. 

Finally,  the  worst  case  scenario  model  gives  the  most  reasonable  results  for  this 
problem.  Both  the  lower  bounds  for  estimates  of  reliabilities  and  95  %  lower  confidence 
bounds  for  these  lower  bounds  are  ordered  according  to  the  total  number  of  failures. 
Thus,  the  results  from  the  worst  case  scenario  were  used  to  implement  the  bonus  system 
for  pyrotechnic  devices. 

After  getting  the  estimate  of  reliability  and  95  %  lower  confidence  bound  for  each 
case,  we  tried  to  design  a  bonus  system  to  improve  the  quality  of  pyrotechnic  devices. 
We  used  95  %  lower  confidence  bounds  instead  of  the  estimated  overall  reliabilities  to 
decide  whether  to  give  bonuses.  Two  sampling  plans  for  giving  the  bonus  to  the  man- 
ufacturers were  considered: 

■  Single  Sampling  Bonus  System 

■  Double  Sampling  Bonus  System. 

We  simulated  two  sampling  bonus  systems  to  see  the  difference  between  them. 
We  concluded  that  the  double  sampling  system  is  more  sensitive  than  single  sampling 
bonus  system. 

To  formulate  an  approximate  bonus  percentage  as  a  function  of  lower  confidence 
bound  for  bonus  (  LCBFB  ),  the  probability  of  passing  manufacturer  test  Pu  and  the 
probability  of  passing  the  environment  tests  P2,  we  used  regression  analysis  with 
GRAFSTAT. 
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Because  the  sampling  plan  results  in  so  much  censoring,  the  only  reasonable  esti- 
mates of  overall  reliability  that  we  obtained  were  actually  lower  bounds.  This  makes  the 
bonus  system  conservative  in  the  sense  that  a  bonus  might  not  be  awarded  what  it  is 
deserved.  Therefore,  if  a  bonus  system  is  to  be  implemented,  a  more  comprehensive 
sampling  plan  needs  to  be  devised  which  allows  estimation  of  R.  A  simple  solution  to 
this  problem  can  be  to  apply  all  environmental  tests  to  the  same  sample  which  would 
give  a  measure  of  dependence  between  environmental  tests. 
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APPENDIX  A.     PROGRAM  MLEA 


PROGRAM  MLEA 


JL. 

* 
* 
* 
* 
* 

* 
* 
* 

* 


THIS  IS  A  FORTRAN  PROGRAM  TO  CALCULATE  THE  RELIABILITY  OF  AN 
ITEM  AFTER  EXPOSURE  TO  SEVERAL  ENVIRONMENTS  WITH   INDEPENDENCE 
ASSUMPTION  WHICH  IS  DESCRIBED  IN  CHAPTER  I  .  THE  PROGRAM  READS 
5000  SUCCESS  VECTORS,  WHICH  ARE  RANDOMLY  GENERATED  BY  PROGRAM 
RANVEC,  FROM  AN  INPUT  FILE  CALLED  SUCVECT  ONE  BY  ONE.  AFTER 
CALCULATION  IT  THEN  PROGRAM  WRITES  ESTIMATED  RELIABILITIES  TO 
AN  OUTPUT  FILE  CALLED  RESULT. 


VARIABLES 

SOM    :  NUMBER  OF  SUCCESSFUL  ITEMS  IN  MANUFACTURER  TEST. 
SOTH   :  NUMBER  OF  SUCCESFUL  ITEMS  IN  TEMPERATURE  AND  HUMIDITY 

TEST. 
SOV    :  NUMBER  OF  SUCCESSFUL  ITEMS  IN  VIBRATION  TEST. 
SOA    :  NUMBER  OF  SUCCESSFUL  ITEMS  IN  ALTITUDE  TEST. 
R1H    :  ESTIMATED  PROBABILITY  OF  PASSING  FROM  MANUFACTURER  TEST. 
R2H    :  ESTIMATED  PROBABILITY  OF  PASSING  FROM  TEMPERATURE  AND 

HUMIDITY  TEST. 
R3H    :  ESTIMATED  PROBABILITY  OF  PASSING  FROM  VIBRATION  TEST 
R4H    :  ESTIMATED  PROBABILITY  OF  PASSING  FROM  ALTITUDE  TEST 
RHMLE  :  ESTIMATED  RELIABILITY  OF  ITEM  AFTER  EXPOSURE  TO 

SEVERAL  ENVIRONMENT  TESTS. 
N     :  SAMPLE  SIZES  FOR  EACH  TEST 
X     :  NUMBER  OF  SUCCESSFUL  ITEMS  IN  EACH  TEST. 
FLAG   :  INDICATOR  VARIABLE  FOR  DETERMINING  EASY  AND  HARD  CASE. 
R1MAX  :  R1H  VALUE  WHICH  MAXIMIZES  LIKELIHOOD. 
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* 


* 
* 


R2H  VALUE  WHICH  MAXIMIZES  LIKELIHOOD. 

R3H  VALUE  WHICH  MAXIMIZES  LIKELIHOOD. 

R4H  VALUE  WHICH  MAXIMIZES  LIKELIHOOD. 

HIGHEST  MAXIMUM  LIKELIHOOD  VALUE. 

LIKELIHOOD  VALUES  AT  THE  END  OF  EACH  HARD  CASE. 

TOTAL  SAMPLE  SIZES  IN  EACH  HARD  CASE. 

REORGANIZED  PROBABILITY  OF  MANUFACTURER  TEST. 

REORGANIZED  PROBABILITY  OF  TEMP.  AND  HUMIDITY  TEST. 

REORGANIZED  PROBABILITY  OF  VIBRATION  TEST. 

REORGANIZED  PROBABILITY  OF  ALTITUDE  TEST. 

JLJLJLJLJL.JLJLJLJLO.J'.JLJL  J-JL  JL  J,  JL   JL  JL  .1.  X  JN  .JL  iJN  ,JL  ,U  ■  JN  .Jt^.  ,U.U.t..U.UmU.U^m^mmU^JmM^.,Ji..Umt^U^».U^.U.l^,t^,Ul1^  JL  JL  JLJL  A  JL  JL  JLJLJL 

TYPE  DECLARATION 

REAL  SOM(5000) ,SOTH(5000) ,SOV(5000) ,SOA(5000) ,R1H,R2H,R3H,R4H, 
+RHMLE(5000),N(4),X(4),R1MAX,R2MAX,R3MAX,R4MAX,LMAX,V(7),Q1H,Q2H, 
+Q3H,Q4H,NTOT,L(7) 

INTEGER  I, J 

LOGICAL  FLAG(4) 


JL 

R2MAX 

JL 

R3MAX 

* 

R4MAX 

JL 

LMAX 

JL 

L 

JL 

V 

JL 

Q1H 

* 

Q2H 

* 

Q3H 

JL 

Q4H 

i-A-ft 


READING  SUCCESS  VECTORS  FROM  SUSVECT  FILE 

DO  60  1=1,5000 

READ(7,*)  SOM(I),SOTH  ( I) ,SOV(I) ,SOA(I) 


•*******■ 


,-A-A-A-AvV  *  ft  A  *  A  A  ft  ft  .V  A  A  -A  A  A  A  A  A  A  A  A  A- .'.  A  J.  A  A  A  ft  A  *  * ft  ft  A  *  A  A  ft  A  ft-A-A-*  ft  ft  ft  A 


FILES  FOR  READING  AND  WRITING 

CALL  EXCMS  ('FILEDEF  7  DISK  SUCVECT  DATA  Al') 

CALL  EXCMS  ('FILEDEF  16  DISK  RESULT  DATA  Al') 

A  ft  A  A  ft  A  A  A  A  A  A  ft  Aftftftftftft  A  A  A  A  A  A  AAA  AAA  AAA  Aftftftftft-A-A-AAAftftftAftAftftftft  A 

INITIALIZATION  OF  SAMPLE  SIZES 

N(l)=  20.0 

N(2)=  20.0 

N(3)=  32.0 

N(4)=  20.0 

INITIALIZATION  OF  FLAG  VARIABLES 

DO  10  J=2,4 

FLAG(J)=  .FALSE. 


JL.JL.JLJUJL.JL 


78 


10  CONTINUE 

*  NUMBER  OF  SUCCESSES  IN  EACH  TEST 
X(l)=  SOM(I) 
X(2)=  SOTH(I) 
X(3)=  SOV(I) 
X(4)=  SOA(I) 

7%  ¥%T\  777S  #»  *»  *»  «»  *»  *»  7w»  »7T;v  **  #»  #v  *%  #»  «i  7V  7*  *»#»**  /*  #»*»**»*#****»***»  *»  #»  **  **  #♦  #»  #*  *»  ?»  «i 

*  CHECK  OPERATION  FOR  EASY  CASE 

IF  ((X(l)/20.00).GE.(X(2)/20.00))  THEN 

FLAG(2)=  .TRUE. 
END  IF 
IF  ((X(l)/20.00).GE.(X(3)/32.00))  THEN 

FLAG(3)=  .TRUE. 
END  IF 
IF  ((X(l)/20.00).GE.(X(4)/20.00))  THEN 

FLAG(4)=  .TRUE. 
END  IF 
IF  (  FLAG(2).AND.FLAG(3).AND.FLAG(4)  )  THEN 

*  *  A  A-  A-  *-A-A-  -A-  A  Me  *  Me  A  A  A  A  *  A  *  A  A  *  A  ft  J<******1rk********ic*it-k*******iric****i^*-k 

*  CALCULATIONS  IN  EASY  CASE 

R1H=  X(l)/20.0 
R2H=  X(2)/20.0 
R3H=  X(3)/32.0 
R4H=  X(4)/20.0 
Q1H=  R1H 
Q2H=  R2H/Q1H 
Q3H=  R3H/Q1H 
Q4H=  R4H/Q1H 

RHMLE(I)=  Q1H*Q2H*Q3H*Q4H 
GO  TO  50 
END  IF 

*  CALCULATIONS  IN  HARD  CASES 

*  INITIALIZATION 
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DO  30  J=l,7 
V(J)=  0.0 
L(J)=  0.0 

30      CONTINUE 
LMAX=   0.0 
NTOT=   0.0 

JLJLJL.JLJLJUJfcO.JU.JUJL.J-JL.^^J^J^JLJI-JLJL-lJL  JLJUJL  JL  JL  JL  JL  jh  JN  ,U.t),,U.tmm1.m1~m1mmUJ.mtmJ~  ^LJUJLi  JL  JL  JL  JLJf^J^  JUJLJUJUJLJL.JLJLJL  JL  JUJL 

•  CASE   1     R1H.LT.1.0     AND     R1H  =  R2H  =  R3H  =  R4H 

DO  40   J=l,4 

V(l)=  V(1)+X(J) 

NTOT=  NTOT+N(J) 
40        CONTINUE 

R1H=  V(l)/NTOT 

L( 1 )=   ( (R1H)**V( 1) )*( (NTOT-V( 1) ) /NTOT)**(NTOT-V( 1) ) 

IF   (L(l).GT.LMAX)    THEN 

LMAX=  L(l) 

R1MAX=  R1H 

R2MAX=  R1H 

R3MAX=  R1H 

R4MAX=  R1H 
END   IF 

ft  Or*****  A  A  A  A  A  A,  A  A  A  A  A  A  A  A  A  A  A  A  A  A- A  ft  -A ft  A-A  A  A -A ******  A  ft  A  *  ft  *  A  ************** 

CASE   2  R1H.LT.1.0,   R1H  =  R2H  =  R3H    ,    R4H   IS  BETWEEN  0.0   AND  R1H 
IF   (   X(4)/20.0.NE.1.0    )   THEN 
V(2)=  X(l)+X(2)+X(3) 
NTOT=  N(l)+N(2)+N(3) 
R1H=  V(2)/NTOT 
IF((X(4)/20.00).LE.R1H)    THEN 

L(2)=   ((R1H)**V(2))*((NT0T-V(2))/NT0T)**(NT0T-V(2))*((X(4) 
+/20.0)**X(4))*(((20.0-X(4))/20.0)**(20.0-X(4))) 
IF   (L(2).GT.LMAX)    THEN 
LMAX=  L(2) 
R1MAX=  R1H 
R2MAX=  R1H 
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R3MAX=  R1H 
R4MAX=  X(4)/20.0 
END   IF 
END   IF 

END   IF 

*  CASE   3   R1H.LT.1.0,    R1H   =   R2H   =   R4H    ,    R3H    IS   BETWEEN   0.0   AND  R1H 
IF   (    X(3)/32.0.NE.1.0    )    THEN 

V(3)=  X(l)+X(2)+X(4) 

NTOT=  N(l)+N(2)+N(4) 

R1H=  V(3)/NTOT 

IF   ((X(3)/32.00).LE.R1H)   THEN 

L(3)=    ((RlH)**V(3))*((NTOT-V(3))/NTOT)**(NTOT-V(3))*((X(3) 
+/32.0)**X(3))*(((32.0-X(3))/32.0)**(32.0-X(3))) 
IF   (L(3).GT.LMAX)   THEN 
LMAX=  L(3) 
R1MAX=  R1H 
R2MAX=  R1H 
R3MAX=  X(3)/32.0 
R4MAX=  R1H 
END   IF 
END   IF 
END   IF 

*-*..*..*    ».i  i  I  i  t  t_t_.t_»  t  ■  I  t  I  ■  I  t.i.t-.t.i.-i..t..t..t..ti.t.tt-«.t«  t . .  t  «  ■  t  «  ■  ■  i  ■ .  »  t  .  t  ■  ■  » .  »  «  .  t  «  i  t..t  t..t..t. 

*  CASE  4  R1H.LT.1.0,    R1H  =  R3H  =  R4H    ,    R2H   IS   BETWEEN   0.0   AND  R1H 
IF   (   X(2)/20.0.NE.1.0    )    THEN 

V(4)=  X(l)+X(3)+X(4) 

NTOT=  N(l)+N(3)+N(4) 

R1H=  V(4)/NT0T 

IF   ((X(2)/20.00).LE.R1H)   THEN 

L(4)=   ((R1H)**V(4))*((NT0T-V(4))/NT0T)**(NT0T-V(4))*((X(2) 
+/20.0)**X(2))*(((20.0-X(2))/20.0)**(20.0-X(2))) 
IF   (L(4).GT.LMAX)   THEN 
LMAX=  L(4) 
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R1MAX=  R1H 
R2MAX=  X(2)/20.0 
R3MAX=  R1H 
R4MAX=  R1H 
END  IF 
END  IF 
END  IF 

*  CASE  5  R1H.LT.1.0,  R1H  =  R2H  R3H  ,  R4H  ARE  BETWEEN  0.0  AND  R1H 

IF  ((X(3)/32.0.NE.1.0).AND.(X(4)/20.0.NE.1.0))   THEN 
V(5)=  X(l)+X(2) 
NTOT=  N(l)+N(2) 
R1H=  V(5)/NTOT 
IF((X(3)/32.00).LE.R1H.AND.(X(4)/20.00).LE.R1H)  THEN 

L(5)=  ((R1H)**V(5))*((NT0T-V(5))/NT0T)**(NT0T-V(5))*((X(3) 
V32.0)**X(3))*((32.0-X(3))/32.0)**(32.0-X(3))*((X(4)/20.0)** 
+x(4))*((20.0-X(4))/20.0)**(20.0-X(4)) 
IF  (L(5).GT.LMAX)  THEN 
LMAX=  L(5) 
R1MAX=  R1H 
R2MAX=  R1H 
R3MAX=  X(3)/32.0 
R4MAX=  X(4)/20.0 
END  IF 
END  IF 
END  IF 

*  CASE  6  R1H.LT.1.0,  R1H  =  R3H  R2H  ,  R4H  ARE  BETWEEN  0.0  AND  R1H 
IF  ((X(2)/20.0.NE.1.0).AND.(X(4)/20.0.NE.1.0))   THEN 

V(6)=  X(l)+X(3) 
NTOT=  N(l)+N(3) 
R1H=  V(6)/NTOT 

IF((X(2)/20.00).LE.R1H.AND.(X(4)/20.00).LE.R1H)  THEN 
L(6)=  ((RlH)**V(6))*((NTOT-V(6))/NTOT)**(NTOT-V(6))*((X(2) 
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+/20.0)**X(2))*((20.0-X(2))/20.0)**(20.0-X(2))*((X(4)/20.0)** 
+x(4))*((20.0-X(4))/20.0)**(20.0-X(4)) 

IF  (L(6).GT.LMAX)  THEN 
LMAX=  L(6) 
R1MAX=  R1H 
R2MAX=  X(2)/20.0 
R3MAX=  R1H 
R4MAX=  X(4)/20.0 
END  IF 
END  IF 
END  IF 

*     CASE  7  R1H.LT.1.0,  R1H  =  R4H   R2H  ,  R3H  ARE  BETWEEN  0.0  AND  R1H 
IF  ((X(2)/20.0.NE.1.0).AND.(X(3)/32.0.NE.1.0))   THEN 
V(7)=  X(l)+X(4) 
NTOT=  N(l)+N(4) 
R1H=  V(7)/NTOT 

IF((X(2)/20.00).LE.R1H.AND.(X(3)/32.00).LE.R1H)  THEN 
L(7)=  ((RlH)**V(7))*((NTOT-V(7))/NTOT)**(NTOT-V(7))*((X(3) 
V32.0)**X(3))*((32.0-X(3))/32.0)**(32.0-X(3))*((X(2)/20.0)**X(2)) 
+*((20.0-X(2))/20.0)**(20.0-X(2)) 

IF  (L(7).GT.LMAX)  THEN 
LMAX=  L(7) 
R1MAX=  R1H 
R2MAX=  X(2)/20.0 
R3MAX=  X(3)/32.0 
R4MAX=  R1H 
END  IF 
END  IF 
END  IF 
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*  CALCULATION  OF  RHMLE  VALUE  WITH  RESPECT  TO  CASE  WHICH  HAS  LARGEST 

*  MAXIMUM  LIKELIHOOD  IN  HARD  CASE 
R1H=  R1MAX 

R2H=  R2MAX 

R3H=  R3MAX 

R4H=  R4MAX 

Q1H=  R1H 

Q2H=  R2H/Q1H 

Q3H=  R3H/Q1H 

Q4H=  R4H/Q1H 

RHMLE(I)=  Q1H*Q2H*Q3H*Q4H 

*  WRITING  AFTER  EACH  CALCULATION 
50   WRITE  (16,*)  RHMLE(I) 

60   CONTINUE 
STOP 
END 
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APPENDIX  B.     PROGRAM  RANVEC 


PROGRAM   RANVEC 


■1. 
* 

* 


* 

a. 


THIS  IS  THE  PROGRAM  TO  GENARETE  RANDOM  NUMBERS  FROM  BINOMIAL 
DISTRIBUTION. THE  PROGRAM  READS  PROBABILITIES  OF  BEING  SUCCESSFUL 
IN  FOUR  TESTS  INTERACTIVELY.  IT  GENERATES  UNIFORMLY  DISTRIBUTED 
RANDOM  NUMBERS  WITH  THESE  PROBABILITIES  ACORDING  TO  SAMPLE  SIZE 
OF  EACH  TEST.  FOR  EACH  TEST  PROGRAM  COUNTS  UNIFORMLY  DISTRIBUTED 
RANDOM  NUMBERS,  WHICH  HAVE  GREATER  THAN  OR  EQUAL  PROBABILITY 
WITH  RESPECT  TO  THE  GIVEN  PROBABILITY  FOR  THAT  TEST. TOTAL  COUNTS 
GIVE  US  SUCCESSFUL  ITEM  NUMBERS  FOR  EACH  TEST.  PROGRAM  UPDATES 
SEEDS  AND  CALLS  SUBROUTINE  RANNUM  IN  EACH  ITERATION.  THE  PROGRAM 
GENERATES  5000  SUCCESS  VECTOR  AND  WRITES  THEM  TO  AN  OUTPUT  FILE 
CALLED  SUSVECT. 


VARIABLES 


PSIM 

PSITH 

PSIV 

PSIA 

NUM1 

NUM2 

NUM3 

NUM4 

X 


PROBABILITY   OF  SUCCESS    IN  MANUFACTURER   TEST. 

PROBABILITY  OF  SUCCESS   IN  TEMP.    AND  HUMIDITY  TEST. 

PROBABILITY  OF  SUCCESS   IN  VIBRATION  TEST 

PROBABILITY  OF  SUCCESS   IN  ALTITUDE  TEST 

COUNTER  FOR  MANUFACTURER  TEST 

COUNTER  FOR  TEMP.    AND  HUMIDITY  TEST 

COUNTER  FOR  VIBRATION  TEST 

COUNTER  FOR  ALTIUTDE  TEST 

BINOMIALLY  DISTRIBUTED  RANDOM  NUMBER  FOR  MANUFACTURER 

TEST 

BINOMIALLY  DISTRIBUTED  RANDOM  NUMBER  FOR  TEMP. AND 

HUMUDITY  TEST 
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* 
* 


* 


c 

D 

ISEED 

KSEED 

LSEED 

MSEED 


BINOMIALLY  DISTRIBUTED  RANDOM  NUMBER  FOR  VIBRATION  TEST 

BINOMIALLY  DISTRIBUTED  RANDOM  NUMBER  FOR  ALTITUDE  TEST 

UNIFORMLY  DISTRIBUTED  RANDOM  NUMBER  FOR  MANUFACTURER 

TEST 

UNIFORMLY  DISTRIBUTED  RANDOM  NUMBER  FOR  TEMP. AND 

HUMIDITY  TEST 

UNIFORMLY  DISTRIBUTED  RANDOM  NUMBER  FOR  VIBRATION  TEST 

UNIFORMLY  DISTRIBUTED  RANDOM  NUMBER  FOR  ALTITUDE  TEST 

SEED  NUMBER  FOR  MANUFACTURER  TEST 

SEED  NUMBER  FOR  TEMP.  AND  HUMIDITY  TEST 

SEED  NUMBER  FOR  VIBRATION  TEST 

SEED  NUMBER  FOR  ALTITUDE  TEST 


TYPE  DECLARATION 

REAL  PSIM,PSITH,PSIV,PSIA,NUM1,NUM2,NUM3,NUM4,X(5000),V(5000) 
+,Y(5000),Z(5000),A,B,C,D 
INTEGER  ISEED, KSEED, LSEED, MSEED 

INITIALIZATION 
ISEED  =  45267 
KSEED  =  113234 
LSEED  =  435 
MSEED  =  1 

*****  A  A  A  A  A  A  A  A  A  A  ft  A  A  *■  A-  ft  *****M"A-******A***A**ti*^^ 

READING  TEST  PROBABILITIES  OF  BEING  SUCCESSFUL  IN  EACH  TEST 

WRITE (*,*)' PLEASE  WRITE  PROBABILITY  OF  BEING  SUCCESSFUL  IN 
+  MANUFACTURER  TEST' 

READ  (*,*)  PSIM 

WRITE(*,*)' PLEASE  WRITE  PROBABILITY  OF  BEING  SUCCESSFUL  IN 
+  TEMPERATURE  AND  HUMIDITY  TEST* 

READ  (*,*)  PSITH 

WRITE(*,*)' PLEASE  WRITE  PROBABILITY  OF  BEING  SUCCESSFUL  IN 
+  VIBRATION  TEST* 

READ  (*,*)  PSIV 
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WRITE(*,*)'PLEASE  WRITE   PROBABILITY  OF  BEING   SUCCESSFUL   IN 
+  ALTITIDE  TEST* 
READ   (*,*)    PSIA 

wi—  -1-  JU  ttg  ~L*  Jm  — I"  »*—  — !■  «J-»  -i—  J«  -*-  ■  I  j  » I  -  -Jj  -t»  -Jj  aJL  -X-  —L.  Jm  — 1»  «J—  —^»  **■»  **-  »■*—  «|jg  ~*—  ■fc  ttg  'fff  »■)■  7^  -*-  *f£*E£  'TPTCfCvCwC4it4&K*ipt$C  «1»  J— -1-«X- J— «JL.«JL.  -J—  «JL.  -J-»  «JL.  -J—  -J—  J-  JL>  -J-  J— 

*  OPENING   AN   OUTPUT  FILE   TO  WRITE   THE   RESULTS 
CALL  EXCMS    ( 'FILEDEF    16   DISK   SUSVECT  DATA   Al') 

*  GENERATION 

DO   50   J=l, 5000,1 

*  INITIALIZATION   IN  EACH   ITERATION 

NUM1=  0.0 
NUM2=  0.0 
NUM3=  0.0 
NUM4=  0.0 

*  SEEDS  UPDATATION   IN  EACH   ITERATION 

ISEED=ISEED+17 
KSEED=KSEED+1356 
LSEED=LSEED+1 
MSEED=MSEED+789 

*  GENERATION  OF  SUCCESSFUL  NUMBER  OF  ITEMS   IN  MANUFACTURER  TEST 

DO   10   1=1,20,1 

CALL  RANNUM   (1,ISEED,0 .0, 1 .0,0 .0,A) 
IF   (   A.LT.PSIM   )   THEN 

NUM1=NUM1+1 
END   IF 
X(J)   =  NUM1 
10  CONTINUE 

*  *A  A  AAA  AAA  A  A  A  A  A  A  A  A Aft  ***JHHHHinllnHH^^ 

*  GENERATION  OF  SUCCESSFUL  NUMBER  OF   ITEMS   IN  TEMP.    AND  HUM.    TEST 

DO  20  K=l,20,l 

CALL  RANNUM   ( 1,KSEED,0.0 , 1. 0,0 . 0 ,B) 
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IF  (  B.LT.PSITH  )  THEN 

NUM2=NUM2+1 
END  IF 
V(J)  =  NUM2 
20        CONTINUE 
ft     ***********■**-. 

*  GENERATION  OF  SUCCESSFUL  NUMBER  OF  ITEMS  IN  VIBRATION  TEST 

DO  30  L=l,32,l 

CALL  RANNUM  ( 1 ,LSEED,0 .0 , 1 . 0 ,0 .0 ,C) 
IF  (  C.LT.PSIV  )  THEN 

NUM3=NUM3+1 
END  IF 
Y(J)  =  NUM3 
30        CONTINUE 

*  GENERATION  OF  SUCCESSFUL  NUMBER  OF  ITEMS  IN  ALTITUDE  TEST 

DO  40  M=l,20,l 

CALL  RANNUM  ( 1 ,MSEED,0 . 0 , 1 . 0 ,0 . 0 ,D) 
IF  (  D.LT.PSIA  )  THEN 

NUM4=NUM4+1 
END  IF 
Z(J)  =  NUM4 
40       CONTINUE 

*  A  A  A  A  A  A  frk  Aft  *rtc*  ft  ft  ft  *  ft  ft  *  ft  *  ft  A  ft  *  ft  *  ft  *  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  A  ft  ft  A  A  ft  A  A  A  A******* 

*  WRITING  THE  RESULTS  TO  AN  OUTPUT  FILE  AS  4  TUPLE 

WRITE  (16,1)  X(J),V(J),Y(J),Z(J) 
1       FORMAT  (1X,F12.7,4X,F12.7,4X,F12.7,4X,F12.7,4X) 

50   CONTINUE 
STOP 
END 

*  A  ft  ft  A  ft  A  A  A  A  A  A  A  ft  ft'  ft  A  ft  ft  A  A  A  A  A  A  ft  A-  A-  A  A-  ft  -A  A  A  ft  ft  ft  A  -A  -A-  A  A  A-  A  ft  A  A  ft  ft  ft  A  A  A  A  A  A  A  A  A  ft-A-ft  ftft-ftft 

SUBROUTINE  RANNUM(DISTN,  SEED,  RPARM1,  RPARM2,  IPARM,  X) 

*  THIS  SUBROUTINE  IS  A  PART  OF  SIMUTIL  FORTRAN  WHICH  IS  WRITTEN 
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*  BY  DR.  M.  P.  BAILEY.   THIS  SUBROUTINE  PROVIDES  AN  INTERFACE  WITH 

*  THE  LLRANDOMII   ROUTINES  PROVIDED  IN  THE  NONIMSL  LIBRARY.   THE 

*  PARAMETER  REQUIRMENTS  AND  CALLING  PROCEDURES  ARE  AS  FOLLOWS: 

*  DISTN  =  DISTRIBUTION  TYPE  YOU  WANT  TO  SELECT  AN  INTEGER  BETWEEN  1 

*  AND  7. 

*  SEED  =  THE  RANDOM  NUMBER  SEED  YOU  WISH  TO  USE. 

*  RPARM1,  RPARM2,  AND  IPARM  ARE  REAL  AND  INTEGER  PARAMETERS  PASSED 

*  TO  THE  ROUTINE  WITH  MEANINGS  WHICH  VARY  WITH  THE  TYPE  OF  DISTRI_ 

*  BUTION  YOU  DESIRE. 

*  X  =  THE  RETURNED  RANDOM  NUMBER,  IT  IS  ALWAYS  REAL. 

*  DISTRIBUTION  NUMBERS  AND  THE  ASSOCIATED  PARM  DEFINITIONS 

*  1- -UNIFORM  ON  THE  INTERVAL  RPARM1  TO  RPARM2 . 

*  2- -NORMAL  WITH  MEAN  RPARM1  AND  VARIANCE  RPARM2. 

*  3- -EXPONENTIAL  WITH  RATE  RPARM1. 

*  4--C0UCHY  WITH  A  =  RPARM1  AND  B  =  RPARM2. 

*  5- -GAMMA  WITH  SHAPE  RPARM2  AND  RATE  RPARM1. 

*  6--POISSON  WITH  RATE  RPARM1. 

*  7--GEOMETRIC  WITH  P  =  RPARM1. 

*  TYPE  DECLARATION 

REAL  RPARM1,RPARM2,X,TEMP,VARIAT(1) 
INTEGER  DISTN,  SEED,  IPARM,  N 

*  *A- A-*A^*A- A  A-^^^  AAA  A  A  A  A  A  Kirk  A  A  A  *  A  A  A  A  A  A  A  A  A  A  A*** 

IF  (DISTN. LE.O. OR. DISTN. GT. 8)  THEN 

WRITE(10,  *)  'ILLEGAL  CALL  TO  RANNUM,  BAD  DISTN' 
STOP 

ENDIF 

GOTO  (10,  20,  30,  40,  50,  60,  70),  DISTN 

*  GENERATE  A  UNIFORM  BETWEEN  RPARM1  AND  RPARM2 
10   CONTINUE 

IF  (RPARM1  -  RPARM2.EQ.0)  THEN 

WRITE(10,  *)  'ILLEGAL  EQUAL  RPARMS  IN  RANNUM' 
STOP 
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ENDIF 

IF   (RPARM1.GT.RPARM2)    THEN 

TEMP  =  RPARM1 

RPARM1   =  RPARM2 

RPARM2   =  TEMP 
ENDIF 

CALL  LRND(SEED,    VARIAT,    1,    1,    0) 

VARIAT(l)    =  RPARM1   +   (RPARM2   -   RPARM1)   *  VARIAT(l) 
GOTO   80 

J-  ifi    1      f      It      I      II      Pii'i    f  ■  ■  W . .  ■  ■  ■  *    il   iTi    I     fn.Hi    I.    *     '■■*ii,.i'    i'i-i*    ■'■■*    i*      '  ■■*■  i*     jj    -*■    [U  JL  Jj  tJU  ■/--*-  -tfy-/t-  /-   '--/    /  -\--'     *■■*    i*    ■*.    *    ■*     *-   *    i'ii'    i'ii'ii*    ■'    .*.i1i 

*  GENERATE   A  NORMAL  WITH  MEAN   RPARM1   AND   STDDEV  RPARM2 
20        CALL  LNORM(SEED,    VARIAT,    1,    1,    0) 

VARIAT(l)    =    (VARIAT(l)    *  RPARM2)    +  RPARM1 
GOTO   80 

JL.  JU.JUJUJUJLJL»JUJUJUJL.JUJ^JUJUJLiJUJU.J^JU 

*  GENERATE  AN  EXPONENTIAL  WITH   RATE    (1/MEAN)    RPARM1 
30        CONTINUE 

IF    (RPARM1.EQ.0)    THEN 

WRITE (10,   *)    'ILLEGAL  ZERO  RATE   IN  RANNUM* 

STOP 
ENDIF 

CALL  LEXPN(SEED,    VARIAT,    1,    1,    0) 
VARIAT(l)   =  VARIAT(l)    /   RPARM1 
GOTO  80 

*  GENERATE  A  COUCHY  WITH  A  =  RPARM1   AND  B   =  RPARM2 
40        CONTINUE 

IF   (RPARM2.LE.0)   THEN 

WRITE(10,    *)    'ILLEGAL  COUCHY   SPREAD   IN  RANNUM,    B  =    ' ,RPARM2 

STOP 
ENDIF 

CALL  LCCHY(SEED,    VARIAT,    1,    1,    0) 
VARIAT(l)   =   (VARIAT(l)   *  RPARM2)    +  RPARM1 
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GOTO  80 


*     GENERATE  GAMMA  WITH  SHAPE  RPARM2  AND  RATE  RPRAM1 
50   CONTINUE 

IF  (RPARM1.LE.0)  THEN 

WRITE(10,  *)  'ILLEGAL  NONPOSITIVE  GAMMA  RATE  IN  RANNUM ' 

STOP 
ENDIF 
IF  (RPARM2.LE.0)  THEN 

WRITE(10,  *)  'ILLEGAL  SHAPE  PARAMETER  IN  RANNUM* 

STOP 
ENDIF 

CALL  LGAMA(SEED,  VARIAT,  1,  1,  0,  RPARM2) 
VARIAT(l)  =  VARIAT(l)  *  (1.0  /  RPARM1) 
GOTO  80 

t     GENERATE  POISSON  WITH  RATE  RPRAM1 
60    CONTINUE 

IF  (RPARM1.LE.0)  THEN 

WRITE (10,  *)  'ILLEGAL  POISSON  RATE  IN  RANNUM' 
STOP 
ENDIF 

CALL  LPOIS(SEED,  VARIAT,  1,  1,  0,  RPARM1) 
GOTO  80 

■  GENERATE  GEOMETRIC  WITH  P  =  RPRAM1 

70    CONTINUE 

IF  (RPARM1.LE.0)  THEN 

WRITE(10,  *)  'ILLEGAL  GEOM  PROB  IN  RANNUM' 
STOP 
ENDIF 

CALL  LGEOM(SEED,  VARIAT,  1,  1,  0,  RPARM1) 
GOTO  80 
80   CONTINUE 
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X  =  VARIAT(l) 
END 
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APPENDIX  C.  PROGRAM  SORT 


* 


PROGRAM  SORT 


THIS  IS  THE  SORTING  PROGRAM.  PROGRAM  USES  BUBBLE  SORT  ALGORITHM. 
PROGRAM  READS  ESTIMATED  RELIABILITIES  FROM  AN  INPUT  FILE  CALLED 
RESULT.  IT  SORTS  FROM   SMALLEST  TO  LARGEST,  AND  WRITES  IN  TO  AN 
OUTPUT  FILE  CALLED  FRESULT  WITH  95  %  LOWER  CONFIDENCE  BOUND. 

VARIABLES 

A      :  ESTIMATED  RELIABILITY 

FLAG   :  INDICATOR  VARIABLE  TAKES  VALUE  '  OF  '  AND  '  OFF  ' 


*  TYPE  DECLARATION 
CHARACTER  FLAG*3 
REAL  A(5000) 
INTEGER  I,N,J 

*  OPENING  AN  INPUT  AND  AN  OUTPUT  FILE 

CALL  EXCMS  ('FILEDEF  9  DISK  RESULT  DATA  Al') 
CALL  EXCMS  ('FILEDEF  15  DISK  FRESULT  DATA  Al') 

*  READING  ESTIMATED  RELIABILITIES 
DO  10  1=1,5000 

READ(9,*)  A(I) 
10   CONTINUE 

*  SORTING  OPERATION 
N=I-1 
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DO  30   I=N,2,-1 
FLAG=*OFF' 
DO   20   J=1,I-1 

IF    (A(J).GT.A(J+1))    THEN 
TEMP=A(J) 
A(J)=A(J+1) 
A(J+1)=TEMP 
FLAG='ON' 
END   IF 
20  CONTINUE 

IF    (FLAG.EQ. 'OFF')    THEN 

GO   TO   40 
END   IF 
30        CONTINUE 
40        CONTINUE 

»  *»  i»  *»  *»  /»  r«  *»  *»  *»  «*  K«C«C«I  R  R  *»'>*»#»  II  it  K  II  II  /»  '»  *»  »»  *»  #%  '*  **  *»  /*  *»  #*  *»  **  #»  #»  **  *»  *»  #»  +%   #»  #*  *»  *»  **  r*  #»  *»  #*  /*  *■%   *-»  *»  t\    r%  1 

fc  WRITING   THE   RESULTS    IN   ASCENDING   ORDER   TO   OUTPUT   FILE 

DO  50    1=1,5000 

WRITE    (15,*)    A(I) 
50        CONTINUE 

WRITE   (15,1)   A(250) 
1  FORMAT    (///,15X,'95   %  LOWER   CONFIDENCE   BOUND    IS ' , 1X.F12 . 7) 

STOP 
END 
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APPENDIX  D.     PROGRAM  INITIAL 


-JL.  -J-_l_-I LJLJ.  J-  J I.XJ.J I 1-~ 


PROGRAM    INITIAL 

*  THIS   PROGRAM,    CALCULATES      INITIAL     GUESSES   FOR   PARAMETERS    IN 

*  LOGLINEAR   MODEL   BY  MEANS   OF   PROGRAM   PARAM   WHICH    IS    IN   APPENDIX   E 

*  PROGRAM   SUPPLIES   PARTIAL   SUMS   OF  EXPECTATION,    TO   SOLVE   EQUATIONS 

*  IN      PROGRAM      PARAM.    IT  READS    INTERACTIVELY   NUMBER   OF   FAILURES    IN 

*  TESTS.    PROGRAM  WRITES   RESULTS   TO  AN   OUTPUT  FILE   CALLED   EXPECT. 

*  VARIABLES 

NUMBER   OF  FAILURES    IN  MANUFACTURER   TEST 

NUMBER   OF  FAILURES    IN  TEMPERATURE   AND   HUMIDITY   TEST 

NUMBER   OF  FAILURES    IN  VIBRATION   TEST 

NUMBER   OF  FAILURES    IN  ALTITUDE   TEST 

SUCCESS   RATIO  FOR  MANUFACTURER   TEST 

SUCCESS   RATIO  FOR  TEMPERATURE   AND   HUMIDITY   TEST 

SUCCESS  RATIO  FOR  VIBRATION  TEST 

SUCCESS  RATIO  FOR  ALTITUDE  TEST 

FAILURE  RATIO  FOR  MANUFACTURER  TEST 

FAILURE  RATIO  FOR  TEMPERATURE  AND  HUMIDITY  TEST 

FAILURE  RATIO  FOR  VIBRATION  TEST 

FAILURE  RATIO  FOR  ALTITUDE  TEST 

PARTIAL  SUM  OF  EXPECTED  NUMBERS   IN  CELLS  FOR  EQUATION    1 

PARTIAL  SUM  OF  EXPECTED  NUMBERS   IN  CELLS  FOR  EQUATION   2 

PARTIAL  SUM  OF  EXPECTED  NUMBERS   IN  CELLS  FOR  EQUATION   3 

PARTIAL  SUM  OF  EXPECTED  NUMBERS   IN  CELLS  FOR  EQUATION  4 


JU 

FOM 

* 

FOTH 

JL 

FOV 

* 

FOA 

* 

PI 

* 

P2 

JU 

P3 

* 

P4 

* 

Ql 

JU 

Q2 

* 

Q3 

JU 

Q4 

XI 

JU 

X2 

JU 

X3 

JU 

X4 

*  TYPE  DECLARATION 
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REAL  F0M,F0TH,F0V,F0A,P1,P2,P3,P4,Q1,Q2,Q3,Q4,PR0D,X1,X2,X3,X4 


*  OPENING  AN  OUTPUT  FILE 

CALL  EXCMS    ('FILEDEF    13   DISK  EXPECT  DATA  Al') 

*  READING   NUMBER   OF  FAILURES    IN   EACH   TEST    INTERACTIVELY 
WRITE(*,*) 'PLEASE   WRITE   NUMBER   OF  FAILURES    IN   MANUFACTURER   TEST* 
READ(*,*)   FOM 

WRITE (*,*)' PLEASE   WRITE   NUMBER   OF  FAILURES    IN  TEMP. AND   HUM. TEST' 

READ(*,*)   FOTH 

WRITE(*,*)' PLEASE   WRITE   NUMBER   OF  FAILURES    IN  VIBRATION   TEST' 

READ(*,*)    FOV 

WRITE(*,*)'PLEASE   WRITE   NUMBER   OF  FAILURES    IN  ALTITUDE   TEST' 

READ(*,*)    FOA 

*  CALACULATION  OF  SUCCESS   RATIOS 
P1=(20.0-FOM)/8.00 
P2=(20.0-FOTH)/8.00 
P3=(32.0-F0V)/8.00 
P4=(20.0-FOA)/8.00 

*  AAAAAA-AA-AAAA"AAAAAAAAAA*AAA^^ 

*  CALCULATION  OF  FAILURE  RATIOS 
Q1=FOM/8.00 

Q2=FOTH/8.00 

Q3=FOV/8.00 

Q4=FOA/8.00 

*- » ■  .  t  -J--J  -  ,1,.  -t , ,  t   I   I   t  *    .  i      t      I   I   i   ■   I   t   t   t._f..t.j..t..i..t.t..t..t.t..t..t   «   ■  .  t . . «  . .  t   t .  t .  i   ti.t.tt   t   t   i   t  _  t   t  .  i  .  t   t.i.t..t_.t.i   ■   t_i._i_ 

*  CALCULATION  PARTIAL  SUMS  OF  EXPECTATION  FOR  PROGRAM  PARAM 
PR0D=2 . 0*P3+2 . 0*Q3+2 . 0*P4+2 . 0*Q4 

Xl=4 . 0*Pl+4 . 0*P2+PROD 
X2=4 . 0*Pl+4 . 0*Q2+PROD 
X3=4 . 0*Q 1+4 . 0*P2+PROD 
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X4=4 . 0*Ql+4 . 0*Q2+PROD 

*     WRITING  THE  RESULTS  TO  AN  OUTPUT  FILE 
WRITE  (13,*)  X1,X2,X3,X4 

STOP 
END 
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APPENDIX  E.     PROGRAM  PARAM 


■k 

* 
* 


PROGRAM   PARAM 


-I 1 1 U-J_J I 1-Jt.U-JL 


■LXJ^XJLUJ^AJ^^AJ^AA^JLAXJLAAAAAAJLA 


j-j-j-u.-i_j_jL.a--i--i--i-j-u--i LJ.XJ.J.J.. 


THIS    IS   THE   PROGRAM   TO   CALCULATE   PARAMETERS   OF  LOGLINEAR   MODEL 
WITH    IMSL   SUBROUTINE.    IT   TAKES   PARTIAL   SUMS  FROM   PROGRAM    INITIAL 
AND   SOLVES      FOUR      NONLINEAR   EQUATIONS,    WHICH   HAVE   FOUR   UNKNOWNS. 
THE   PROGRAM  USES   AN    IMSL   SUBROUTINE   CALLED  DNEQNF  TO   SOLVE   THIS 
EQUATION.IT  WRITES   SOLUTIONS   OF  EQUATIONS   TO  OUTPUT  FILE   CALLED 
PARAM  DATA. 

VARIABLES 

ITMAX  :    MAXIMUM    ITERATION  NUMBER. 

N  :    PARAMETER 

XGUESS  :    INITIAL  GUESS  FOR  FOUR  NONLINEAR  EQUATIONS. 

F  :    NONLINEAR  EQUATIONS. 

TYPE  DECLARATION 

PARAMETER   (N=4) 
REAL*8   ERRREL 
INTEGER   ITMAX.N 


* 
* 


INTEGER  K 

REAL*8     FNORM,X(N),XGUESS(N) 

EXTERNAL  ACN 

OPENING  A  FILE  FOR  WRITING  RESULTS 


-JUJUJL  -JL.  J-  -l—l.^l—l-JUJL. LJLJL  JUJ— JL 
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CALL  EXCMS  ('FILEDEF  9  DISK  PARAM  DATA  A') 

*  INITIAL  GUESS 
DATA  XGUESS/3 . 5D0 , 3 . 5D0 , 3 . 5D0 , 3 . 5D0/ 

*  INITIALIZATION 
ERRREL  =  0.0001D0 
ITMAX  =  10000 

*  CALLING  OF  IMSL  SUBROUTINE 
CALL  DNEQNF  (  ACN, ERRREL, N, ITMAX ,XGUESS ,X ,FNORM) 

*  RESULTS 

WRITE  (9,1)  (X(K),K=l,N),FNORM 
1   FORMAT('THE  SOLUTION  TO  THE  SYSTEM  IS ',/, 'X=( ' ,4F8 . 2 ,')',/, 'WITH 
+FNORM=' ,F8.2,//) 
END 

SUBROUTINE  ACN  (X,F,N) 

*  »***********  A  A  A  A  ****  ****iriTi!izirk**iiitiTifit^itiziTicif*ifis*iiiT^if**ifitisizi;it**iT-k* 

*  VARIABLES 

*  X   :  INITIAL  GUESS 

*  F   :  NONLINEAR  EQUATIONS 

*  TYPE  DECLARATION 
REAL*8  X(N),F(N) 
INTEGER  N 

*  1  ST  EQUATION 
F(1)=DEXP(X(1))*(DEXP(X(2)+X(3)+X(4))+ 

+  (1/DEXP(X(2)+X(3)+X(4)))+ 

+  (DEXP(X(2)+X(3))/DEXP(X(4)))+ 

+  (DEXP(X(2)+X(4))/DEXP(X(3)))+(DEXP(X(3) 

+  +X(4))/DEXP(X(2)))+(DEXP(X(2))/DEXP(X(3)+X(4)))+ 
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+     (DEXP(X(3))/DEXP(X(2) 

+      +X(3)))+(DEXP(X(4))/DEXP(X(2)+X(3))))-LOG(32.00) 

2  ST  EQUATION 
F(2)=DEXP(X(1)+X(2))*(DEXP(X(3)+X(4))+ 

+     (1/DEXP(X(3)+X(4)))+(DEXP(X(3))/ 

+     DEXP(X(4)))+(DEXP(X(4))/DEXP(X(3))))-LOG(23.00) 

3  ST  EQUATION 
F(3)=DEXP(X(1)+X(3))*(DEXP(X(2)+X(4))+ 

+     (1/DEXP(X(2)+X(4)))+(DEXP(X(2))/ 

+     DEXP(X(4)))+(DEXP(X(4))/DEXP(X(2))))-LOG(23.00) 

4  ST  EQUATION 
F(4)=DEXP(X(1)+X(4))*(DEXP(X(2)+X(3))+ 

+      (1/DEXP(X(2)+X(3)))+(DEXP(X(2))/ 

+      DEXP(X(3)))+(DEXP(X(3))/DEXP(X(2))))-LOG(14.00) 

RETURN 

END 
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APPENDIX  F.     PROGRAM  LLMDEP 


__L.J — t_-t_J — U-l__) — L.  J-  J — l_-J — I_J 1 I-  J.a_  J : 


PROGRAM   LLMDEP 


* 
* 


* 


THIS    IS   THE      FORTRAN   PROGRAM   TO   CALCULATE   THE  RELIABILITY 
OF      PYROTECHNIC   DEVICE.       IT     ASSUMES      THAT  THERE    IS   A   DEPENDENCE 
BETWEEN   MANUFACTURER  AND  ENVIRONMENT  TESTS.    EXPECTATION-MAXIMIZA. 
TION   ALGORITHM    IS   USED    IN   THIS   MODULE.         ALGORITHM      STARTS   WITH 
INITIAL  GUESSES   FOR   PARAMETERS      AND  ESTIMATES   EXPECTATIONS.    IT 
RECALCULATES  CELL   PROBABILITIES   AND  UPDATES   EXPECTATIONS   UNTIL 
IT   CONVERGES.      AN    ITERATIVE   NEWTON   AND  RAPHSON   PROCEDURE    IS   USED 
DURING   UPDATATION   OF  CELL   PROBABILITIES.      THIS   PROCEDURE    IS   DONE 
BY   A   SUBROUTINE  NAMED  UCPROB . 


VARIABLES 


■^^MU^UJUJUJUJUJUJUJU^J 


***• 


-l.J-.r-«i*-I,.^l-.J,-t.^.^l..l..l..l..t..J..J..I.^...t..T..I.JI-J-^.J.-,l— t.J^JL.JU-JUJL,U^U.U— UJ— L 


J--L.JL.J 


JU 

FOM 

■k 

FOTH 

* 

FOV 

* 

FOA 

* 

MU 

JL 

LP1 

* 

LP2 

JL 

LP3 

* 

LP4 

a. 

TETHA 

* 

RHMLE 

* 

MPPPP 

* 

MPOPO 

* 

NUMBER  OF  FAILURES    IN  MANUFACTURER  TEST. 
NUMBER  OF  FAILURES   IN  TEMP.    AND   HUMIDITY  TEST. 
NUMBER  OF  FAILURES   IN  VIBRATION  TEST. 
NUMBER  OF  FAILURES   IN  ALTITUDE   TEST. 
OVERALL  MEAN. 

MEAN  EFFECT  OF  MANUFACTURER  TEST. 
MEAN  EFFECT  OF  MANUFACTURER  TEST. 
MEAN  EFFECT  OF  MANUFACTURER  TEST. 
MEAN  EFFECT  OF  MANUFACTURER  TEST. 
TWO  WAY  INTERACTION  TERMS 
RELIABILITY  OF  DEVICE 

CELL  FREQUENCY  WITH  RESPECT  TO  TESTS  RESULTS. 
SUM  OF  CELL  FREQUENCIES   WHICH   HAVE   PASSED  DEVICES 
MANUFACTURER  AND  VIBRATION  TEST. 
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* 

MPOOO 

* 

* 

IP 

* 

FP 

* 

Y 

* 

A,B,C 

* 

EXPPPP 

* 

* 

* 

FLAG 

Vr 

RHMLE 

SUM  OF  CELL  FREQUENCIES  WHICH  HAVE   PASSED  DEVICES 

MANUFACTURER  TEST. 

INITIAL  PROBABILITY  VECTOR 

UPDATED    (FINAL)    PROBABILITY   VECTOR 

CELL  EXPECTATION   VECTOR 

SOME   TERMS   TO  MAKE   THE   CALCULATIONS   EASY. 

EXPECTED  NUMBER   OF  DEVICES    IN   CELL  WHICH 

HAS   A   RESULTANT  VECTOR    (    P   P   P   P    )    IN   MANUFACTURER, 

TEMPERATURE,      VIBRATION  AND  ALTITUDE  TEST  RESPECTIVELY 

INDICATOR   VARIABLE   OF  CONVERGENCE  FOR   PARAMETERS. 

ESTIMATED  RELIABILITY   OF   PYROTECHNIC   DEVICE. 

*-JL.-A-JL.-X-»U»L.  J.O.,LJ.  JL  JL  JL  J- ,^-l.JUa-JLJU  JU-l.Jt-JU  JL.U.U ■!■■!■  ■!■  ,J,.I,,I,.J— 1- 0,0. -UP.  JLO— I- -l-JU.l-a..l- .1,  J.  0,  0,0.0,  JUO-JU-l-JU-l—JU-J—l—l,,!..!, 
ft  ft  ft  ft  ft  ft  «»  t*    i\    r%    it  ft  *»  li  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  f%  ft  ft  ft  ft  ft  ft  f»  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  ft  f%  ft  ft  ft  ft  ft  ft 

*  TYPE   DECLERATION 

PARAMETER    (K=10000) 

LOGICAL  FLAG (16) 

INTEGER    I 

REAL  FOM,    FOTH,    FOV,    FOA ,    LP1,    LP2 ,    LP3 ,    LP4 ,    TETHA 

REAL  MPPPP,MPPPF,MPPFP,MPPFF, 
+  MPFPP,MPFPF,MPFFP,MPFFF, 

+  MFPPP,MFPPF,MFPFP,MFPFF, 

+  MFFPP,MFFPF,MFFFP,MFFFF 

REAL  EXPPPP(K) ,EXPPPF(K) ,EXPPFP(K) ,EXPPFF(K) , 
+  EXPFPP(K) ,EXPFPF(K) ,EXPFFP(K) ,EXPFFF(K) , 

+  EXFPPP(K) ,EXFPPF(K) ,EXFPFP(K) ,EXFPFF(K) , 

+  EXFFPP(K) ,EXFFPF(K) ,EXFFFP(K) ,EXFFFF(K) 

REAL  MPOOO,    MFOOO,    MPPOO,    MPOPO ,    MPOOP,    A,    B,    C 

REAL    IP(16,1),FP(16,1),X(16) 

REAL  Y(16) 

COMMON    /    PROB    /   Y 

*  -  •*■  -'•  '-  ~  "  ~  -  •'•  -  -  -  ~  -  "  *•  •'■  ~  *  ~  i'.  >'-  -  •'• >'-  •'■  •'•  ■'■  J.AAA  .'. ;'.  A  *  *  A  .',  A  J. )'.  -A  A-  *  A  ,', A  A ,',  A  A  ;.  a  ,', ,',  )', ,',  -i,.  a  a  a ,'. 

CALL  EXCMS   ('FILEDEF   15   DISK  END  DATA  Al    '    ) 

*  il  'I  .'.  ■■'.  L-  .'.  1   1  1  ,'.  ■!   r'.   ,',  ,?,    I,  ;',   ,i    '    .'   ■'   \   I    I    <    I    r    l    l    '   ji    i    i    i    i   ^   i^  ^  ;_  ^  ^   i%  ^  ^  ^   .^  ^  ^  ^   ,%  |  |  j>  ^   .^   .    ^  m   .^   ._  m  »  >  |   . 

*  INITIALIZATION 

RHMLE   =0.0 
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EPS  =  0.001 
DO  10  L=l,16 

FLAG (L)=. FALSE. 
10   CONTINUE 

*  READING  THE  NUMBER  OF  FAILURES  IN  EACH  TEST  INTERACTIVELY 
WRITE(*,*) 'PLEASE  ENTER  THE  #  OF  FAILURES  IN  MANUFACTURER  TEST' 
READ(*,»)  FOM 

WRITE(*,*) 'PLEASE  ENTER  THE  #  OF  FAILURES  IN  TEMP.  AND  HUM.  TEST' 

READ(*,*)  FOTH 

WRITE(*,*) 'PLEASE  ENTER  THE  #  OF  FAILURES  IN  VIBRATION  TEST' 

READ(*,*)  FOV 

WRITE(*,*) 'PLEASE  ENTER  THE  #  OF  FAILURES  IN  ALTITUDE  TEST' 

READ(*,*)  FOA 

*  READING  THE  INITIAL  GUESS  FOR  EACH  CELL  IN  HYPOTHETICAL 

*  CONTINGENCY  TABLE 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MPPPP' 

READ(*,*)  MPPPP 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MPPPF* 

READ(*,*)  MPPPF 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MPPFP' 

READ(*,*)  MPPFP 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MPPFF' 

READ(*,*)  MPPFF 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MPFPP' 

READ(*,*)  MPFPP 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MPFPF' 

READ(*,*)  MPFPF 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MPFFP' 

READ(*,*)  MPFFP 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MPFFF* 

READ(*,*)  MPFFF 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MFPPP' 
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READ(*,*)  MFPPP 

WRITE  (*,*)'PLEASE  ENTER  INITIAL  GUESS  FOR  MFPPF' 

READ(*,*)  MFPPF 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFPFP ' 

READ(*,*)  MFPFP 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFPFF' 

READ(*,*)  MFPFF 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFFPP' 

READC*,*)  MFFPP 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFFPF' 

READ(*,*)  MFFPF 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFFFP' 

READC*,*)  MFFFP 

WRITE  (*,*) 'PLEASE  ENTER  INITIAL  GUESS  FOR  MFFFF' 

READC*,*)  MFFFF 

*  CALCULATION  OF  INITIAL  PROBABALITIES  USING  CELL  FREQUENCIES 
IP(1,1)  =  MPPPP/92.00 
IP(2,1)  =  MPPPF/92.00 
IP(3,1)  =  MPPFP/92.00 
IP(4,1)  =  MPPFF/92.00 
IP(5,1)  =  MPFPP/92.00 
IP(6,1)  =  MPFPF/92.00 
IP(7,1)  =  MPFFP/92.00 
IP(8,1)  =  MPFFF/92.00 
IP(9,1)  =  MFPPP/92.00 
IPC  10,1)=  MFPPF/92.00 
IP(11,1)=  MFPFP/92.00 
IP(12,1)=  MFPFF/92.00 
IP(13,1)=  MFFPP/92.00 
IP( 14,1)=  MFFPF/92.00 
IPC  15,1)=  MFFFP/92.00 
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IP(16,1)=  MFFFF/92.00 

*  DETERMINATION  OF   INITIAL  FREQUENCIES  FOR  LIKELIHOOD  ESTIMATION 
Y(l)    =  MPPPP 

Y(2)  =  MPPPF 
Y(3)  =  MPPFP 
Y(4)  =  MPPFF 
Y(5)  =  MPFPP 
Y(6)  =  MPFPF 
Y(7)  =  MPFFP 
Y(8)  =  MPFFF 
Y(9)  =  MFPPP 
Y(10)=  MFPPF 
Y(ll)=  MFPFP 
Y(12)=  MFPFF 
Y(13)=  MFFPP 
Y(14)=  MFFPF 
Y(15)=  MFFFP 
Y(16)=  MFFFF 

*  CALL  A   SUBROUTINE  WHICH  UPDATES  CELL  PROBABILITIES  USING  NEWTON 

*  AND  RAPHSON  PROCEDURE 
CALL  UCPROB(    IP,FP   ) 

*-  1--I  -  .  t  _-!__■   1   t   ■   t   1..I..1..1..1   1.1.  .1..1..I..1.  1  .  1  .  .*..*.     1.1.  1.1. .1.1   1.1   1.1.  ■  »-  .  t  .  1  .  .1..1..1..1.  .1.1..*   1   1  .  1..  t.  .  t  ..  1  .  .  1   1..1..1..I..1   1   1   I   1   I   1   I  -t- 

*  UPDATATION  OF  CELL  FREQUENCIES 
MPPPP=FP(1,1)*92.00 
MPPPF=FP(2,1)*92.00 
MPPFP=FP(3,1)*92.00 
MPPFF=FP(4,1)*92.00 
MPFPP=FP(5,1)*92.00 
MPFPF=FP(6,1)*92.00 
MPFFP=FP(7,1)*92.00 
MPFFF=FP(8,1)*92.00 
MFPPP=FP(9,1)*92.00 
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MFPPF=FP(10,1)*92.00 
MFPFP=FP(11,1)*92.00 
MFPFF=FP(12,1)*92.00 
MFFPP=FP(13,1)*92.00 
MFFPF=FP( 14,1)^92.00 
MFFFP=FP(15,1)*92.00 
MFFFF=FP(16,1)*92.00 


■  '■  ■1..fl..f..l..'..l..t..'..'.  .'■■'■  ■'  t  ?  t  t  t  t.  t  .  t  t  t  .  t.  t  t  i  t  t  t  t  t  ^  i  |LJUJU^UJUJUJUJUJUJUJUJUJ|  UJUJUJI  '  '%  t  *  '  '  '  '  '  I  jLJ  ' 


INITIAL  EXPECTATIONS 


EXPPPP(1 
EXPPPF(1 
EXPPFP(1 
EXPPFF( 1 
EXPFPPU 
EXPFPF(1 
EXPFFP(1 
EXPFFF(1 
EXFPPP( 1 

EXFPPFd 
EXFPFP(1 
EXFPFF(1 
EXFFPP(1 
EXFFPF(1 
EXFFFP(1 
EXFFFF( 1 


)=MPPPP 
)=MPPPF 
)=MPPFP 
)=MPPFF 
)=MPFPP 
)=MPFPF 
)=MPFFP 
)=MPFFF 
)=MFPPP 
)=MFPPF 
)=MFPFP 
)=MFPFF 
)=MFFPP 
)=MFFPF 
)=MFFFP 
)=MFFFF 


A  A  A- A  A  A  A  A  A  A  A  A  A  A  -A-A- A  ft  trick***  A  A  A  A  A  A  A  A  A  A  A  A  A  A  ********  ****************** 

MPOOO=MPPPP+MPPPF+MPPFP+MPPFF+MPFPP+MPFPF+MPFFP+MPFFF 

MFO  0  0  =MFPPP+MFPPF+MFPFP+MFPFF+MFFPP+MFFPF+MFFFP+MFFFF 

MPPOO=MPPPP+MPPPF+MPPFP+MPPFF 

MPOPO=MPPPP+MPPPF+MPFPP+MPFPF 

MPOOP=MPPPP+MPPFP+MPFPP+MPFFP 

A=92.00-MPP00 

B=92.00-MP0P0 
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C=92.00-MP00P 

NEXT  EXPECTATIONS 
DO  20  1=2, K 

EXPPPP(I)=(20.00-FOM)*(MPPPP/MP000)+(20.00-FOTH)*(MPPPP/MPP00)+ 

-  (32.00-FOV)*(MPPPP/MP0P0)+(20.00-FOA)*(MPPPP/MP00P) 
IF(ABS(EXPPPP(I)-EXPPPP(I-1)) .LE.EPS)  THEN 

FLAG(1)=.TRUE. 
END  IF 
EXPPPF(I)=(20.00-FOM)*(MPPPF/MP000)+(20.00-FOTH)*(MPPPF/MPP00)+ 

-  (32.00-FOV)*(MPPPF/MP0P0)+FOA*(MPPPF/C) 
IF( ABS(EXPPPF( I ) -EXPPPFC I- 1 ) ) . LE . EPS)  THEN 

FLAG(2)=.TRUE. 
END  IF 

EXPPFP(I)=(20.00-FOM)*(MPPFP/MP000)+(20.00-FOTH)*(MPPFP/MPP00)+ 
■  FOV*(MPPFP/B)+(20.00-FOA)*(MPPFP/MP00P) 
IF( ABS(EXPPFP( I ) -EXPPFP( I- 1 ) ) . LE . EPS)  THEN 

FLAG(3)=.TRUE. 
END  IF 

EXPPFF(I)=(20.00-FOM)*(MPPFF/MP000)+(20.00-FOTH)*(MPPFF/MPP00)+ 
FOV*( MPPFF/B ) +FOA*( MPPFF/C ) 
IF(ABS(EXPPFF(I)-EXPPFF(I-1)).LE.EPS)  THEN 

FLAG(4)=.TRUE. 
END  IF 

EXPFPP(I)=(20.00-FOM)*(MPFPP/MP000)+FOTH*(MPFPP/A)+(32.00-FOV)* 
(MPFPP/MP0P0)+(20.00-FOA)*(MPFPP/MP00P) 
IF(ABS(EXPFPP(I)-EXPFPP( I- 1)). LE.EPS)  THEN 

FLAG(5)=.TRUE. 
END  IF 

EXPFPF(I)=(20.00-FOM)*(MPFPF/MP000)+FOTH*(MPFPF/A)+(32.00-FOV)* 
( MPFPF/MPOPO ) +FOA* ( MPFPF/C ) 
IF(ABS(EXPFPF(I)-EXPFPF(I-1)). LE.EPS)  THEN 

FLAG(6)=.TRUE. 
END  IF 
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EXPFFP(I)=(20.00-FOM)*(MPFFP/MP000)+FOTH*(MPFFP/A)+FOV*(MPFFP/B)+ 
+   (20.00-FOA)*(MPFFP/MP00P) 

IF(ABS(EXPFFP(I)-EXPFFP(I-1)).LE.EPS)  THEN 

FLAG(7)=.TRUE. 
END  IF 

EXPFFF(I)=(20.00-FOM)*(MPFFF/MP000)+FOTH*(MPFFF/A)+FOV*(MPFFF/B)+ 
+  FOA*(MPFFF/C) 

IF(ABS(EXPFFF(I)-EXPFFF(I-1)).LE.EPS)  THEN 

FLAG(8)=.TRUE. 
END  IF 

EXFPPP(I)=FOM*(MFPPP/MF000)+FOTH*(MFPPP/A)+FOV*(MFPPP/B)+FOA* 
+   (MFPPP/C) 

IF(ABS(EXFPPP(I)-EXFPPP(I-1)).LE.EPS)  THEN 

FLAG(9)=.TRUE. 
END  IF 

EXFPPF(I)=FOM*(MFPPF/MF000)+FOTH*(MFPPF/A)+FOV*(MFPPF/B)+FOA* 
+   (MFPPF/C) 

IF(ABS(EXFPPF(I)-EXFPPF(I-1)).LE.EPS)  THEN 

FLAG(10)=.TRUE. 
END  IF 

EXFPFP(I)=FOM*(MFPFP/MF000)+FOTH*(MFPFP/A)+FOV*(MFPFP/B)+FOA* 
+   (MFPFP/C) 

IF(ABS(EXFPFP(I)-EXFPFP(I-1)).LE.EPS)  THEN 

FLAG(11)=.TRUE. 
END  IF 

EXFPFF ( I ) =FOM* ( MFPFF/ MFO  0  0 ) +FOTH* ( MFPFF/ A ) +FOV* ( MFPFF/ B ) +FOA* 
+   (MFPFF/C) 

IF(ABS(EXFPFF(I)-EXFPFF(I-1)).LE.EPS)  THEN 

FLAG (12)=. TRUE. 
END   IF 

EXFFPP(I)=FOM*(MFFPP/MFOOO)+FOTH*(MFFPP/A)+FOV*(MFFPP/B)+FOA* 
+      (MFFPP/C) 

IF(ABS(EXFFPP(I)-EXFFPP(I-1)).LE.EPS)    THEN 
FLAG(13)=.TRUE. 
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END  IF 

EXFFPF ( I ) =FOM* ( MFFPF/ MFO  0  0 ) +FOTH* ( MFFPF/ A ) +FOV* ( MFFPF/ B ) +FOA * 
+   (MFFPF/C) 

IF(ABS(EXFFPF(I)-EXFFPF(I-1)).LE.EPS)  THEN 
FLAG(14)=.TRUE. 

END  IF 

EXFFFP(I)=FOM*(MFFFP/MFOOO)+FOTH*(MFFFP/A)+FOV*(MFFFP/B)+FOA* 
+   (MFFFP/C) 

IF(ABS(EXFFFP(I)-EXFFFP(I-1)).LE.EPS)  THEN 
FLAG U5)=. TRUE. 

END  IF 

EXFFFF ( I ) =FOM* ( MFFFF/ MFO  0  0 ) +FOTH* ( MFFFF/ A ) +FO V* ( MFFFF/ B ) +FO A* 
+   (MFFFF/C) 

IF(ABS(EXFFFF(I)-EXFFFF(I-1)).LE.EPS)  THEN 
FLAG(16)=.TRUE. 

END  IF 


MPPPP 
MPPPF 
MPPFP 
MPPFF 
MPFPP 
MPFPF 
MPFFP 
MPFFF 
MFPPP 
MFPPF 
MFPFP 
MFPFF 
MFFPP 
MFFPF 
MFFFP 


EXPPPP(I 
EXPPPF(I 
EXPPFP(I 
EXPPFF(I 
EXPFPPd 
EXPFPF(I 
EXPFFP(I 
EXPFFF(I 
EXFPPP(I 
EXFPPF(I 
EXFPFP(I 
EXFPFF(I 
EXFFPP(I 
EXFFPF ( I 
EXFFFP(I 
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MFFFF  =  EXFFFF(I) 

*  CHECK  FOR  THE   STOPING  CONDITION 

IF(   FLAG(1).AND.FLAG(2).AND.FLAG(3).AND.FLAG(4).AND.FLAG(5) 
+      . AND . FLAG ( 6 ) . AND . FLAG ( 7 ) . AND . FLAG ( 8 ) . AND . FLAG ( 9 ) . AND . FLAG (10) 
+      .AND.FLAG(11).AND.FLAG(12).AND.FLAG(13).AND.FLAG(14).AND. 
+     FLAG(15).AND.FLAG(16))   THEN 

*  CALCULATION  OF  FINAL  EXPPPP    (STOPING   CONDITION    IS   MET) 

RHMLE=EXPPPP(I)/92.00 
GO  TO   30 
END   IF 

*  STOPING   CONDITION   IS   NOT  MET.      PROBABILITIES   FOR   THE   NEXT 
NEWTON  AND  RAPHSON   PROCEDURE 

IP(1,1)  =  MPPPP/92.00 
IP(2,1)  =  MPPPF/92.00 
IP(3,1)  =  MPPFP/92.00 
IP(4,1)  =  MPPFF/92.00 
IP(5,1)  =  MPFPP/92.00 
IP(6,1)  =  MPFPF/92.00 
IP(7,1)  =  MPFFP/92.00 
IP(8,1)  =  MPFFF/92.00 
IP(9,1)  =  MFPPP/92.00 
IP( 10,1)=  MFPPF/92.00 
IP(11,1)=  MFPFP/92.00 
IP(12,1)=  MFPFF/92.00 
IP(13,1)=  MFFPP/92.00 
IP( 14,1)=  MFFPF/92.00 
IPC  15,1)=  MFFFP/92.00 
IPC16,1)=  MFFFF/92.00 

*A  A  A  A  A  A  A  A  A  A  A  A  A  ;':  A  A  A  A  A  A  A  A  A  A  ..'  i?  HH    f    I    I    '  J    '■  «  ^■■^■*   !    I    '    r..t    t.  j.  .i„  i    i    i    i    i   i   ■    i  ,i,  i,  i,,i.  i    i    i    i    t 

CALL  UCPROB(IP,FP) 

"   "   "   ''   r»   *w»   .»   *»  rt  ?»  #*  *»  «w«  <t  ^»  /«  *w«  r*  #»  *%  #*  /*  *»  #»  /»  /»  *»   „   rv  *W%  *»   R  R   R   #■*   >»   .»   ,»   ,,   ..   .,   /»  V,   . \  V.   #»   #\   R   .\   ,  .   A   ff  'ii  VV  V%   #\   *  .   ME  «V  V* 
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MPPPP=FP(1,1)*92.00 

MPPPF=FP(2,1)*92.00 

MPPFP=FP(3,1)*92.00 

MPPFF=FP(4,1)*92.00 

MPFPP=FP(5,1)"92.00 

MPFPF=FP(6,1)*92.00 

MPFFP=FP(7,1)*92.00 

MPFFF=FP( 8, 1)*92.00 

MFPPP=FP(9,1)*92.00 

MFPPF=FP(10,1)*92.00 

MFPFP=FP(11,1)*92.00 

MFPFF=FP(12,1)*92.00 

MFFPP=FP(13,1)*92.00 

MFFPF=FP(14,1)*92.00 

MFFFP=FP(15,1)*92.00 

MFFFF=FP(16,1)*92.00 

*  **********^**A"A-A"A-*c**Vc*'^^ 

MPOOO=MPPPP+MPPPF+MPPFP+MPPFF+MPFPP+MPFPF+MPFFP+MPFFF 

MFOOO=MFPPP+MFPPF+MFPFP+MFPFF+MFFPP+MFFPF+MFFFP+MFFFF 

MPPO  0=MPPPP+MPPPF+MPPFP+MPPFF 

MPOPO=MPPPP+MPPPF+MPFPP+MPFPF 

MPO  0  P=MPPPP+MPPFP+MPFPP+MPFFP 

A=92.00-MPP00 

B=92.00-MP0P0 

C=92.00-MP00P 

TrwvrTrTrTrTrTr7r7rTr7rvrTrTrTV7tTrTr 

20      CONTINUE 

*  **fr*  A  A  A  A  A  A  A  A  A  A  A  JHHHHlrtHr^HlnllriHI^^ 

30      WRITE(15,40)FOM,FOTH,FOV,FOA,RHMLE 

40      FORMATC/.SX/CASE'^X.FS^^X^S^^X.FS^^X.FS^./.ISX, 
+*MLE   =    ',F12.7) 

STOP 
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END 


XJ-J.J.J-J-- LJLJ.U.J.J.J-J-J-J-. 


■  J-J — I — l_J-J_-i-J_U_U-Jt--J — UUUUUUULJLAJULJLjL. l.JUJt-0.. LA«UA«LAJUJJULJLJLJULJL 


SUBROUTINE  UCPROB    (IP,FP) 

THIS   SUBROUTINE   UPDATES   CELL   PROBABILITIES   USING  NEWTON   AND 
RAPHSON   PROCEDURE   WHICH    IS   DESCRIBED    IN   SAS    . 


VARIABLES 


-1 1-  J- J 1 1 1 UU.U.-1-. 


'■■>!■■% 


-1-  J- J 1 LJ 1 T L.X..J 1 t.JUJUJUU.JUJL.^L 


***■ 


-J-JL.J. 


J.JL,-! 


r* 


JL 

FOM 

* 

FOTH 

* 

FOV 

JL. 

FOA 

* 

MU 

Jb 

LP1 

* 

LP2 

* 

LP3 

* 

LP4 

* 

TETHA 

* 

RHMLE 

* 

MPPPP 

* 

MPOPO 

* 

* 

MP000 

* 

* 

IP 

* 

FP 

* 

Y 

* 

A,B,C 

* 

EXPPPP 

* 

* 

* 

SIGN 

* 

RHMLE 

NUMBER  OF  FAILURES    IN   MANUFACTURER   TEST. 
NUMBER  OF  FAILURES    IN   TEMP.    AND   HUMIDITY   TEST. 
NUMBER   OF  FAILURES    IN  VIBRATION  TEST. 
NUMBER  OF  FAILURES   IN  ALTITUDE  TEST. 
OVERALL  MEAN. 

MEAN  EFFECT  OF  MANUFACTURER  TEST. 
MEAN  EFFECT  OF  MANUFACTURER  TEST. 
MEAN  EFFECT  OF  MANUFACTURER  TEST. 
MEAN  EFFECT  OF  MANUFACTURER  TEST. 
TWO  WAY  INTERACTION  TERMS 
RELAIBILITY  OF  DEVICE 

CELL  FREQUENCY  WITH  RESPECT  TO  TESTS  RESULTS. 
SUM   OF  CELL  FREQUENCIES  WHICH   HAVE   PASSED  DEVICES 
MANUFACTURER  AND  VIBRATION  TEST. 

SUM  OF  CELL  FREQUENCIES  WHICH  HAVE   PASSED  DEVICES 
MANUFACTURER  TEST. 
INITIAL  PROBABILITY  VECTOR 
UPDATED   (FINAL)    PROBABILITY  VECTOR 
CELL  EXPECTATION  VECTOR 

SOME  TERMS  TO  MAKE  THE  CALCULATIONS  EASY. 
EXPECTED  NUMBER  OF  DEVICES   IN  CELL  WHICH 
HAS  A  RESULTANT  VECTOR   (    P  P  P  P   )    IN  MANUFACTURER, 
TEMPERATURE,      VIBRATION  AND  ALTITUDE  TEST  RESPECTIVELY 
INDICATOR  VARIABLE  OF  CONVERGENCE  FOR   PARAMETERS. 
ESTIMATED  RELIABILITY  OF  PYROTECHNIC  DEVICE. 
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TYPE  DECLERATION 

LOGICAL  SIGN(5) 

REAL   IP(16,1),FP(16,1) 

REAL  F0(15,1),F1(15,1),S(15,15),X(15,5),B0(5,1), 
+B1(5,1),PI0(15,1),PI1(15,1),C(5,5),G(5,1),SINV(15,15),F(15,1), 
+PR|(15,1),PR2(5,1),PR3(15,5),PR4(5,5),SUM,LAST,PI01(16,1),LHE, 
+LHEMAX,RIP(15,1),PR5(15,5),CINV(5,5),DIF(15,1),PR6(5,1),LAMBDA, 

+PI11(16,1),EPS,BLAST(5,1),PLAST(15,1),BNEW(5,1),FLAST(15,1) 
+,U,XT(5,15),PR4INV(5,5) 

REAL   Y(16) 

COMMON   /   PROB   /   Y 

INTEGER    I,J,K 

INITIALIZATION 

DO  50    1=1,5 

SIGN(I)=    .FALSE. 
50      CONTINUE 

DO   70    1=1,5 
DO   60   J=l,5 
S(I,J)=0.0 
60  CONTINUE 

70      CONTINUE 

***********  A  A  A  A  .V  A  A  A**********************  A ******************** 

READING  THE  DESIGN  MATRIX 

CALL  EXCMS   ('FILEDEF  9  DISK  DESIGN   INPUT     Al') 

DO  80   1=1,15 

READ(9,*)   X(I,1),X(I,2),X(I,3),X(I,4),X(I,5) 
80      CONTINUE 
REWIND  9 

.i,f.ir.f.-i-,i..i..i,,t,j,.f.,i,.t..t  .,«  .  t._t_  .  t   f  .  _i   t   i  ._  t  --I-  .t..i-.i  _-  t--i,--,i .  ,t ..  _t  -- 1  -  ,t  -_  t  .  ..r .  t  .  t  .  i    _  t  _t  . .  t  ..  i__  t .  -»  . .  t  .  .  > . .  t  .  .  t  .  t  .  t  .    i.i.    t  .  t  .  t  .  . «  .  .  t . .  i  . .  t  .  .  t  .  i  ..  i .  .  i  .  - 1  . 

INVERSE  OF  VARIANCE  AND  COVERIANCE  MATRIX  FOR   INITIAL  BO 
DO   100    1=1,15 
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DO  90  J=l,15 

IF(I.EQ.J)  THEN 

SINV(I,I)=(IP(I,1)-(IP(I,1)**2.0))*92.00 
END  IF 
IF(I.NE.J)  THEN 

SINV(I,J)=(-IP(J,1))*IP(I,1)*92.00 
END  IF 
90     CONTINUE 
100   CONTINUE 

*  LOGIT  RESPONSE  FUNCTIONS 
DO  110  1=1,15 

F(I,1)  =ALOG(IP(I,l)/IP(16,l)) 
110   CONTINUE 

*  ***********  -A-  *  *  -k  A  rt  *  *  A  ■/>■  A  aWc  A  A-  *  v'v  vV  ■;'.■  ***********  A  A  A  A  ******************** 

*  THE  TRANSPOZE  OF  THE  DESIGN  MATRIX 
CALL  TRNRR  ( 15,5 ,X, 15,5 , 15 ,XT,5) 

*  *AAAA-A-A-AAAAAAAAAAAAAAAAAAA-A-ft-A--A&A-A-A*A-AAAAAAAAAAAAAAA^^ 

*  MATRIX  MULTIPLICATION  PR1=(SINV*F) 

CALL  MRRRR  ( 15 , 15 ,SINV, 15 , 15 , 1 ,F, 15 , 15 , 1.PR1 , 15) 

*  &■**  A  A**-*- AA"A-/rAJ"A^^ 

*  MATRIX  MULTIPLICATION  PR2=(XT*PR1) 

CALL  MRRRR  (5, 15 ,XT,5 , 15 , 1.PR1, 15 ,5 , 1 ,PR2,5) 

*  ****A-**AvVvVA^^ 

*  MATRIX  MULTIPLICATION  PR3=(SINV*X) 

CALL  MRRRR  ( 15 , 15 ,SINV, 15, 15,5 ,X, 15 , 15 ,5 ,PR3, 15) 

*  MATRIX  MULTIPLICATION  PR4=(XT*PR3) 

CALL  MRRRR  (5 , 15 ,XT,5 , 15 ,5 ,PR3, 15 ,5 ,5 ,PR4,5) 

*  AA  A  A  A  A  A  A  A  A  A- A  A  AAA  A  A  A  A  A  A  A  A  frMetctrk  A A-A-A ****************  A  A  A  A  A  KX******* 

*  INVERSE  OF  THE  MATRIX  MULTIPLICATION  PR4=PR4INV 
CALL  LINRG  (5,PR4,5,PR4INV,5) 

*  ******  A"A  A  AAA  A  A  AAA  A  AA^^ 
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*  INITIAL  ESTIMATION  OF  PARAMETERS  BO 

*  MATRIX  MULTIPLICATION  B0=(PR4INV*PR2) 

CALL  MRRRR  (5 ,5 ,PR4INV,5 ,5 , 1 ,PR2 ,5 ,5 , 1 ,BO ,5) 

*  F0=  X  *  BO 

CALL  MRRRR  ( 15 ,5 ,X , 15 ,5 , 1 ,B0 ,5 , 15 , 1 ,F0 , 15) 

*  INITIAL  PROBABILITIES  PI0=(EXP(F0) ) 
DO  120  1=1,15 

PI0(I,1)=EXP(F0(I,1)) 
120   CONTINUE 
SUM=0.0 
DO  130  1=1,15 

SUM=SUM+PI0(I,1) 
130   CONTINUE 

LAST=1.0/(1.0+SUM) 

*  PROBABILITY  MATRIX  WHICH  INCLUDES  16  VALUES   PI01 

*  FOR  THE  INITIAL  ESTIMATE  OF  LIKELIHOOD  ESTIMATION 
DO  140  1=1,15 

PI01(I,1)=PI0(I,1)*LAST 
140   CONTINUE 

PI01(16,1)=LAST 

*  INITIAL  LIKELIHOOD  FOR  NEXT  ITERATION  AT  STEP  BO 
LHE=0 . 0 

DO  150  1=1,16 

LHE=LHE+Y(I)*ALOG(PI01(I,1)) 
150   CONTINUE 

LHEMAX=LHE 

*  REORGANIZED  INITIAL  PROBABILITIES  FOR  UPDATATION  RIP 
DO  160  1=1,15 

RIP(I,1)=IP(I,1) 
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160   CONTINUE 


*  FIRST  ITERATION  IN  NEWTON  AND  RAPHSON  METHOD 
DO  170  1=1,15 

DO  180  J=l,15 

IF(I.EQ.J)  THEN 

SINV(I,I)=(PI0(I,1)-(PI0(I,1)**2.0))*92.00 
END  IF 
IF(I.NE.J)  THEN 

SINV(I,J)=(-PI0(J,l))*PIO(I,l)*92.00 

END  IF 
170     CONTINUE 
180   CONTINUE 

*  ************?W.**A**AAAAAA^^ 

*  MATRIX  MULTIPLICATION  PR5=(SINV*X) 

CALL  MRRRR  (  15  ,  15  ,SINV,  15  ,  15  ,5  ,X  ,  15  ,  15  ,5  ,PR5  ,  15) 

*  ************************************^ 

*  MATRIX  MULTIPLICATION  C=(XT*PR5) 

CALL  MRRRR  (5 , 15 ,XT,5 , 15 ,5 ,PR5 , 15 ,5 ,5 ,C,5) 

*  *********  A  A  A  A  A  A  A  /.  A  A  A  A  A  A  A  A  A  A  AAA  A  A  A  A  A  AA  A**-A-**"A*-h"A***************** 

*  INVERSE  OF  THE  MATRIX  C=CINV 
CALL  LINRG  (5 ,C ,5 ,CINV,5) 

DO  190  1=1,15 

DIF(I,1)=  92.00*(RIP(I,1)-PI0(I,1)) 
190   CONTINUE 

*  **********  A  A  A  A  A  A  A  A  A  A  A  A  A  A  *****************  k-k-k-klrk-fcjtftjc-?;  ■!<&■/*■/:&&&■&&■&•&•& 

*  MATRIX  MULTIPLICATION  FOR  G=(XT*DIF) 

CALL  MRRRR  (5 , 15, XT, 5, 15 , 1 ,DIF, 15 ,5, 1,G,5) 

*  ****************************^ 

*  MATRIX  MULTIPLICATION  (PR6=  CINV*G) 

CALL  MRRRR  (5,5,CINV,5,5 , 1 ,G,5,5 , 1 ,PR6,5  ) 

*  *  ****  A  -A' A*  A  A  A  A  A  A  A  A  A  A  A  A  A  -A  A  ********  A  A  A  *•  ***iricicicic±^/{i"h"A"\-itiririrk*icicicic*ic* 

LAMBDA=1.0 
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200      DO  210    1=1,5 

PR6(I,1)=PR6(I,1)*LAMBDA 
210      CONTINUE 


"  INITIAL   VALUE   FOR   Bl 

DO   220    1=1,5 

B1(I,1)=B0(I,1)-PR6(I,1) 
220      CONTINUE 


*  MATRIX   MULTIPLICATION      F1=(X*B1) 

CALL  MRRRR    ( 15 ,5 ,X, 15 ,5 , 1 ,B1 ,5 , 15 , 1 ,F1 , 15) 


UJUJLJUJUJUJLJJUJU.JUJUJ, 


*  CALCULATION   OF   PROBABILITIES   FOR   Bl   PI1=EFP(F1) 
DO   230    1=1,15 

PI1(I,1)=EXP(F1(I,1)) 
230      CONTINUE 

*  CALCULATION   OF  THE    16   TH   PROBABILITY  VALUE 
SUM=0.0 
DO  240   1=1,15 

SUM=SUM+PI 1(1,1) 
240      CONTINUE 

LAST=1.0/(1.0+SUM) 

JU  .t-.t-.i  ..t-.t.-t.  t._t..t.  t  .t  .»   *  J   t.  *..t..i.  *<»*«««   it   «ttttttttt« 

DO  250   1=1,15 

PI11(I,1)=PI1(I,1)*LAST 
250      CONTINUE 

PI11(16,1)=LAST 

*  INITIAL  LIKELIHOOD  ESTIMATION 
LHE=0.0 

DO  260    1=1,16 

LHE=LHE+Y( I )*ALOG( PI 11(1,1)) 
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260   CONTINUE 


K=0 

EPS=0.001 

IF(LHE.LT.LHEMAX)  THEN 
K=K+1 

IF(K.GT.IO)  THEN 
DO  270  1=1,5 

BLAST(I,1)=B0(I,1) 
270         CONTINUE 

DO  280  1=1,15 

PLAST(I,1)=PI01(I,1) 
280         CONTINUE 
GO  TO  490 
END  IF 

LAMBDA=LAMBDA/2.0 
GO  TO  200 
END  IF 
LHEMAX=LHE 

*  *******  A  A  A  ft  A  A  A  A  A  A  A  A  A  A  A  A  A  A  A  A  ft  A  A  A  A  rt  *  A  A  A  A A  A  J.  A  A  fc******  Ifk+rk-k-k*******-* 
DO  290  1=1,5 

IF(ABS(B0(I,1)-B1(I,1)).LE.EPS)  THEN 

SIGN(I)=.TRUE. 
END  IF 
290   CONTINUE 

*  CHECKING  CRITERIAS 
IF(SIGN(1).AND.SIGN(2).AND.SIGN(3).AND.SIGN(4).AND.SIGN(5))  THEN 

DO  300  1=1,5 

BLAST(I,1)=B1(I,1) 
300      CONTINUE 

DO  310  1=1,15 

PLAST( I, 1)=PI 11(1,1) 
310      CONTINUE 
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GO   TO  490 
END    IF 

*  CRITERIAS  ARE  NOT  MET  THEN  NEW   ITERATIONS 
DO   320    1=1,5 

BLAST(I,1)=B1(I,D 
320      CONTINUE 

DO   330    1=1,15 

PLAST(I,1)=PI1(I,1) 
330      CONTINUE 

340      DO   360    1=1,15 

DO   350   J=l,15 

IF(I.EQ.J)    THEN 

SINV(I,I)=(PLAST(I,1)-(PLAST(I,1)**2.0))    *92.00 
END   IF 
IF(I.NE.J)   THEN 

SINV(I,J)=(-PLAST(J,1))*PLAST(I,1)*92.00 
END   IF 
350  CONTINUE 

360     CONTINUE 

**  **  *»  #v  #w»7w\  #*  T»*»75*»7»  7»  *\  71  «»  if  It  71  7%  7C  7E  IS  H  7f  7»  7*  7*  7%  7*  7f  75  7i  7%  7»  *»  *v  if  7*  *»  7»  7*  7»  **  *»  **7V  #»  <»  **  #v  **  #»  #*  *»  #»  #*  <»  <»  #*  **  #» 

*  MATRIX  MULTIPLICATION  PR3=(SINV*X) 

CALL  MRRRR   ( 15 , 15 ,SINV, 15 , 15 ,5 ,X, 15 , 15 ,5 ,PR3, 15) 

*  MATRIX  MULTIPLICATION  C=(XT*PR3) 

CALL  MRRRR    (5 , 15 ,XT,5 , 15 ,5 ,PR3, 15 ,5 ,5 ,C,5) 

CALL  LINRG   (    5,C,5,CINV,5    ) 

DO  370   1=1,15 

DIF(I,1)=  92.00*(RIP(I,1)-PLAST(I,1)) 


119 


370   CONTINUE 

JU  JLJUJL.JUJLJUJLJL,., 

*  MATRIX  MULTIPLICATION  G=(XT*DIF) 
CALL  MRRRR  (5 , 15 ,XT,5 , 15 , 1,DIF, 15 ,5 , 1,G,5) 

*  MATRIX  MULTIPLICATION  PR6=(CINV*G) 
CALL  MRRRR  (5 ,5 ,CINV,5 ,5 , 1 ,G,5 ,5 , 1 ,PR6,5) 

J-  JL._L.JL.JL.  J_JL.U_-lrJ L,. 

LAMBDA=1.0 

JU  JLJL.JL.JLJL,. JL.-t-JL.JL— JL- 

*  NEW  PARAMETER  ESTIMATES 
380   DO  390  1=1,5 

BNEW( I , 1 )=BLAST( 1,1) - (LAMBDA*PR6 (1,1)) 
390   CONTINUE 

JL  JUJVJUJUJJL^LJUJUJJUJL-JL-J^ 

*  MATRIX  MULTIPLICATION  F1=(X*BNEW) 

CALL  MRRRR  ( 15 ,5 ,X, 15 ,5 , 1 ,BNEW,5 , 15 , 1 ,F1 , 15) 
DO  400  1=1,15 

PI1(I,1)=EXP(F1(I,1)) 
400   CONTINUE 

SUM=0 . 0 

DO  410  1=1,15 

SUM=SUM+PI 1(1,1) 
410   CONTINUE 

LAST=1.0/(1.0+SUM) 

*  1ririrtrJrk-kirirkiricM"k  A  aHe  A  A  A  A1  *  *  *  A-  ft  ft  ft  '!s*ic*1rkirtciTititiriiit1:ir&-ki;iti?i;-ki;-k-kir-k*-k-k-ki>-k-k* 

DO  420  1=1,15 

PI11(I,1)=PI1(I,1)*LAST 
420   CONTINUE 

PI11(16,1)=LAST 

*  LIKELIHOOD  ESTIMATION 
LHE=0.0 
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DO  430    1=1,16 

LHE=LHE+Y(I)*ALOG(PI 11(1,1)) 

430      CONTINUE 


K=0 

EPS=0.001 

IF(LHE.LT.LHEMAX)    THEN 

K=K+1 

IF(K.GT.IO)    THEN 
GO   TO  490 

END   IF 

LAMBDA=LAMBDA/2.0 

GO  TO   380 
END   IF 
LHEMAX   =LHE 

DO  440    1=1,5 

IF(ABS(B1(I,1)-BLAST(I,1)).LE.EPS)   THEN 

SIGN(I)=.TRUE. 
END   IF 
440      CONTINUE 

-L.  -t--i.-y..t-y-.t..i     *.   t   .t ..«. ,!,.«,.,  x    »,  ■«     *.,'  ■  i»  ,».«■■*■■«.   *■■«,.*..'■  .1.  .1.  .1.  .t,  .1   JL  J,  .1.  .I,    I.  J   .1.  .t.if  ..B  i.f     lull  I?  .f   i*.  *nf  ii '■  UJUJUJ«JUJ*JUJUJUJUJUJL 

IF(SIGN(1).AND.SIGN(2).AND.SIGN(3).AND.SIGN(4).AND.SIGN(5))   THEN 
DO  450   1=1,5 

BLAST(I,1)=BNEW(I,1) 
450  CONTINUE 

DO  460   1=1,15 

PLAST(I,1)=PI11(I,1) 
460  CONTINUE 

GO  TO  490 
END   IF 

*  CRITERIAS  ARE  NOT  MET  THEN  NEW   ITERATIONS 

DO  470   1=1,5 
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BLAST(I,1)=B1(I,1) 
470      CONTINUE 

DO  480    1=1,15 

PLAST(I,1)=PI11(I,1) 

480      CONTINUE 
GO  TO   340 


JU  -T..,tr.,t  ,tr-,tr-.>r  i ,,  t  j  j  i  i  t  t  t  t      i  -t  .-t.  t  i..«   i._i„i„t.-t„i_.i  ■■t.-i..t.-t^-_^-.r..^.-r.-^t--^..t-.t.-t--i-^i-^t-^t--t,.t--i--i--t-.t..,t„t.  ,t,  .,*..t..f  ,.!..■ 

490      CONTINUE 

*  MATRIX  MULTIPLICATION   FLAST=(X*BLAST) 
CALL  MRRRR    ( 15 ,5 ,X, 15 ,5 , 1 .BLAST ,5 , 15 , 1 ,FLAST, 15) 

*  CALCULATION  OF  FINAL  PROBABILITIES 
DO   500    1=1,15 

PLAST( I , 1)=EXP(FLAST( 1,1)) 
500      CONTINUE 

*  ******  *******  ***  AAA,  *  **  ************************** 

SUM=0.0 

DO  510    1=1,15 

SUM=SUM+PLAST(I,1) 
510      CONTINUE 

*  A  A  A  A  A  A  A  A  A  A  A  A  A  A  A A  A  A  A  ********************************************* 

LAST=1.0/(1.0+SUM) 
DO  520    1=1,15 

FP( I , 1)=PLAST( I , 1 )*LAST 
520      CONTINUE 

FP(16,1)=LAST 

RETURN 

END 
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APPENDIX  G.     PROGRAM  MLEB 


.JL..I..I.U-JU-U.I..JL.JU>JJU-J 


PROGRAM   MLEB 

J  Aj  A*  Ji!  iT '  *J  *  'JhP  *J 

THIS    IS   THE   PROGRAM   TO  CALCULATE   RELIABILITY   OF  THE   DEVICE 

*  WITH   DEPENDENCE   ASSUMPTION . PROGRAM   ASSUMES   THAT  FAILED   ITEM   FROM 

*  ANY   OF  ENVIRONMENT  TESTS   FAILS      FROM   MANUFACTURER   TEST   TOO.    THIS 

*  IS   WORST  CASE   SCENARIO.    IT  READS   NUMBER   OF   SUCCESSFUL   ITEMS   FROM 

*  AN    INPUT  DATA   CALLED   SUSVECT.    FINALLY   THE   PROGRAM  WRITES   RESULTS 

*  TO  AN   OUTPUT  FILE   CALLED  RESULT. 


* 

* 
* 


VARIABLES 

Ai'ii'i  ' 

SOM  :    NUMBER   OF   SUCCESSFUL   ITEMS    IN  MANUFACTURER   TEST 

SOTH         :    NUMBER   OF   SUCCESSFUL   ITEMS    IN   TEMP.    AND   HUMIDITY   TEST 

SOV  :    NUMBER   OF   SUCCESSFUL   ITEMS    IN  VIBRATION   TEST 

SOA  :    NUMBER   OF   SUCCESSFUL   ITEMS   IN  ALTITUDE   TEST 

R1H  :    ESTIMATED   PROBABILITY  OF   PASSING  FROM  MANUFACTURER  TEST 

R2H  :    ESTIMATED  PROBABILITY  OF  PASSING  FROM  TEMPRATURE  AND 

HUMIDITY  TEST 
R3H  :    ESTIMATED  PROBABILITY  OF  PASSING  FROM  VIBRATION  TEST 

R4H  :   ESTIMATED  PROBABILITY  OF  PASSING  FROM  ALTITUDE  TEST 

RHMLE      :    ESTIMATED  RELIABILITY  OF   ITEM  AFTER  EXPOSURE  TO 

SEVERAL  ENVIRONMENT  TESTS. 
X  :    DUMMY  VARIABLE 

TYPE  DECLARATION 

REAL  SOM(5000) ,SOTH(5000) ,SOV(5000) ,SOA(5000) ,X(4) ,R1H,R2H,R3H, 
+  R4H,RHMLE(5000) 
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INTEGER  I 

*  FILES  FOR  READING  AND  WRITING 

CALL  EXCMS  ( 'FILEDEF  9  DISK  SUSVECT  DATA  Al') 
CALL  EXCMS  ('FILEDEF  17  DISK  RESULT  DATA  Al1) 

*  READING  NUMBER  OF  SUCCESS  IN  EACH  TEST 
DO  10  1=1,5000 

READ(9,*)  SOM(I),SOTH(I),SOV(I) ,SOA(I) 

JU  ■1..f..l..l..r..l..l,.l..l..l..f..t..r..,..*..l..,..1..j,..'..'..f..,).  ■'..'.  X  JL  JL  JL  ,«^  JL  JL  JUA  ^  JL^JLJUAJUJLiAJUJ.JLJL  JUJL  ^  .'.  .*..f..*..l..l..!  ..I..I..I ..!.,»  .'..!. 

*  NUMBER  OF  SUCCESSES  IN  EACH  TEST 

X(l)=  SOM(I)+SOTH(I)+SOV(I)+SOA(I) 

X(2)=  SOTH(I) 

X(3)=  SOV(I) 

X(4)=  SOA(I) 

*  CALCULATIONS   SUCCESS   PROBABILITIES   IN  EACH  TEST 

R1H=  X(l)/184.0 

R2H=  X(2)/40.0 

R3H=  X(3)/64.0 

R4H=  X(4)/40.0 

RHMLE(I)=  R1H  *  R2H  *  R3H  *  R4H 

*  ft  ftft  ft  ft  *  -A  ft  ft  ft  ft  ft  ft  A  ft  A  A  A  A  A'  ft  ft  ft  ft  ft  ft  ft  A  A  A  A  A  A  A  A  *  ft  A  ft  ft  A  A  A  ft  ******************** 

*  WRITING  RESULTS  TO  AN  OUTPUT  FILE  CALLED  RESULT 

WRITE  (17,*)  RHMLE(I) 
10    CONTINUE 
STOP 

END 
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APPENDIX  H.     PROGRAM  BONUS 


PROGRAM   BONUS 


*  THIS  IS  THE  PROGRAM  TO  CALCULATE  BONUS  PERCENTAGE  OF  ANY  FIRM 

*  WHOSE  LONG  RUN  SUCCESS   PROBABILITIES  ARE  KNOWN.  IN  THIS  PROGRAM 

*  IT  USES  2000  SUCCESS  VECTORS,  WHICH  ARE  GENERATED  BY  RANVEC  IN 

*  APPENDIX   B.  THEY  ARE  GENERATED  BY  KNOWN  LONG  RUN  PROBABILITIES 

*  THE  PROGRAM  USES  TWO  DATA  SETS. THEY  ARE  PRECALCULATED  LCB ' S  SETS 

*  .FIRST  DATA  REPRESENTS  FIRST  INSPECTION,  SECOND  DATA  REPRESENTS 

*  SECOND  INSPECTION  LOWER  CONFIDENCE  BOUNDS. 

*  SUCCESS  VECTORS  REPRESENT  OFFERED  LOTS. IT  HAS  A  DETERMINISTIC 

*  BONUS  LINE. PROGRAM  CALCULATES  LCB  OF  OFFERED  LOT, WITH  FIRST  DATA 

*  AND  COMPARES  IT  WITH  LCB  OF  BONUS  LINE. IF  FIRM  LCB  IS  GRATER  THAN 

*  FIRM  GETS  BONUS.  OTHERWISE  FIRM  HAS  A  CHANCE  TO  ONE  MORE  TRY.  IN 

*  SECOND  TRY,  PROGRAM  CUMULATES  SUCCESS  VECTORS  AND  IT  USES  SECOND 

*  DATA  TO  FIND  OUT  LCB  OF  CUMULATED  LOT.  AFTER  THIS  CALCULATION  IT 

*  COMPARES  AGAIN.  FINALLY  IT  COUNTS  NUMBER  OF  TIMES  THAT  THE  FIRM 

*  GETS  THE  BONUS  IN  1000  REPLICATIONS  AND  ESTIMATES  BONUS  PERCENT. 

*  IT  WRITES  RESULTS  TO  AN  OUTPUT  FILE  CALLED  BONUS  DATA. 

*  *****************************^ 

*  VARIABLES 

*  ********  ,V  A  A  A  A  A  A  A  A  A  A A **  A  A  A  A  A  A  ************************************ 

NUMBER  OF  SUCCESS  IN  MANUFACTURER  TEST. 

NUMBER  OF  SUCCESS  IN  TEMPERATURE  AND  HUMIDITY  TEST. 

NUMBER  OF  SUCCESS  IN  VIBRATION  TEST. 

NUMBER  OF  SUCCESS  IN  ALTITUDE  TEST. 

LOWER  CONFIDENCE  BOUND  OF  BONUS  LINE 

LOWER  CONFIDENCE  BOUND 

LCB  VALUES  ARRAY  IN  FIRST  INSPECTION 


*  SOM 

*  SOTH 

*  SOV 

*  SOA 

*  BLINE 

*  LCB 


*     FILCB 
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SELCB 

A,B,C,D 

BFI 

BSI 

BTOT 

PRCT 

COUNT 

SIGN 

FLAG 


LCB  VALUES  ARRAY  IN  SECOND  INSPECTION 

DIMENSIONS  FOR  USE  OF  LCB  DATAS 

NUMBER  OF  TIMES  THAT  FIRM  GAT  BONUS  AFTER  1  ST  INSP. 

NUMBER  OF  TIMES  THAT  FIRM  GAT  BONUS  AFTER  2  ST  INSP. 

TOTAL  NUMBER  OF  TIMES  THAT  FIRM  GAT  BONUS. 

BONUS  PERCENT. 

COUNTER  FOR  1000  REPLICATIONS. 

INDICATOR  OF  ACCEPTANCE  FOR  FIRST  INSPECTION 

INDICATOR  OF  ACCEPTANCE  FOR  SECOND  INSPECTION 


*  TYPE  DECLARATION 
LOGICAL  SIGN(4),FLAG(4) 

REAL  SOM(2000),SOTH(2000),SOV(2000),SOA(2000),BLINE,LCB,BFI,BSI, 
+  BTOT,FILCB(2,2,3,2),SELCB(3,3,5,3),A,B,C,D,PRCT,COUNT 
INTEGER  I,J,K,L 

*  OPENING  FILES  FOR  READING  AND  WRITING 

CALL  EXCMS  (  '  FILEDEF  7  DISK  SUCVECT  DATA  Al'  ) 

CALL  EXCMS  (  '  FILEDEF  8  DISK  FIRST  DATA  Al'  ) 

CALL  EXCMS  (  '  FILEDEF  9  DISK  SECOND  DATA  Al'  ) 

CALL  EXCMS  (  *  FILEDEF  15  DISK  BONUS  DATA  Al'  ) 

*  INITIALIZATION 
COUNT  =1.0 
BFI  =  0.0 
BSI  =  0.0 
BLINE  =  0.9250000 
DO  10  1=1,4 

SIGN(I)=.TRUE, 
FLAG(I)=.TRUE, 
10    CONTINUE 

JU  JU  JL  JU  JU  ■  jU  JU  JU  ■  jU  .^1  ■  J  ■  ■!-  JL  JU  JL  JL  JL  J 

READING  FIRST  INSPECTION  LOWER  CONFIDENCE  BOUNDS  FROM  DATA  FILE 
DO  50  1=1,2 
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DO  40  J=l,2 

DO  30   K=l,3 

DO   20   L=l,2 

READ(8,*)   FILCB(I,J,K,L) 

20 

CONTINUE 

30 

CONTINUE 

40 

CONTINUE 

50 

CONTINUE 

*  READING   SECOND    INSPECTION   LOWER   CONFIDENCE   BOUNDS  FROM   DATA   FILE 
DO   90    1=1,3 

DO   80   J=l,3 

DO   70   K=l,5 

DO   60   L=l,3 

READ(9,*)    SELCB(I,J,K,L) 
60  CONTINUE 

70  CONTINUE 

80  CONTINUE 

90        CONTINUE 

#*  n  #v  **  *\  4<*  *%  t\  '» *\  »v  7v7»  n  #*  44  iii  vi  *v  n  mni  n  n  n  ttttt*  *»  7w»  **  *»  /» **  **  #»  *»  *»  **  #»7*  <*  *%  *»?»#*  <»  *»7w»tw»  **  n  #»?»  *»**#*  *»  #* 

*  READING   SUCCESS   PROBABILITY  OF  FIRM    IN  EACH   TEST 

WRITE(*,*)' WRITE   THE   PROBABILITY   OF   SUCCESS    IN  MANUFACTURER   TEST* 
READ   (*,*)    PSIM 

WRITE(*,*)'WRITE  THE   PROBABILITY  OF  SUCCESS   IN  TEMPERATURE  AND 
+HUMIDITY  TEST' 
READ   (*,*)   PSITH 

WRITE(*,*)'WRITE  THE  PROBABILITY  OF  SUCCESS   IN  VIBRATION  TEST' 
READ   (*,*)   PSIV 

WRITE(*,*)'WRITE  THE   PROBABILITY  OF  SUCCESS   IN  ALTITUDE  TEST' 
READ   (*,*)   PSIA 

*  READING     NUMBER  OF  SUCCESFUL   ITEMS 

*  FOR  FIRST   INSPECTION  AND  SECOND   INSPECTION 
DO   100   1=1,2000 
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READ(7,*)  SOM(I),SOTH(I),SOV(I),SOA(I) 
100   CONTINUE 

*  INSPECTIONS  BEGIN 
DO  120  N=l,2000,2 

*  SUBSCRIPT  DEFINITION   FOR  FIRST  INSPECTION 

A=INT(20.0-SOM(N))+1 
B=INT(20.0-SOTH(N))+1 
C=INT(32.0-SOV(N))+1 
D=INT(20.0-SOA(N))+1 

*  CHECK  FOR  ACCEPTANCE  OF  FIRST  OFFERED  LOT 

IF  (A.GT.2)  THEN 

SIGN(1)=. FALSE. 
END  IF 
IF  (B.GT.2)  THEN 

SIGN(2)=. FALSE. 
END  IF 
IF  (C.GT.3)  THEN 

SIGN(3)=.FALSE. 
END  IF 
IF  (D.GT.2)  THEN 

SIGN(4)=. FALSE. 
END  IF 
IF( .NOT. (SIGN(l) .AND.SIGN(2) .AND.SIGN(3) .AND.SIGN(4) ))  THEN 

SIGN(1)=.TRUE. 

SIGN(2)=.TRUE. 

SIGN(3)=.TRUE. 

SIGN(4)=.TRUE. 

GO  TO   112 
END  IF 

*  *falHlriHHIr:)HlriHi^ 

*  AFTER  FIRST  INSPECTION  DETERMINATION  OF  LCB  OF  FIRM 
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LCB=FILCB(A,B,C,D) 


IS  LOT  WORTH  WHILE  FOR  GETTING  BONUS  ? 
IF  (  LCB.GT.BLINE  )  THEN 

BFI=BFI+1.0 

GO  TO   110 
END  IF 


*  SUBSCRIPT  DETERMINATION  FOR  SECOND  INSPECTION 

A=A+INT(20.0-SOM(N+1)) 
B=B+INT(20.0-SOTH(N+1)) 
C=C+INT(32.0-SOV(N+1)) 
D=D+INT( 20 . 0-SOA(N+1) ) 

*  CHECKING   RESULTS  OF  SECOND  TEST  SERIES  ABOUT  ACCEPTANCE 

IF  (A.GT.3)  THEN 

FLAG(1)=. FALSE. 
END  IF 
IF  (B.GT.3)  THEN 

FLAG(2)=.FALSE. 
END  IF 
IF  (C.GT.5)  THEN 

FLAG(3)=.FALSE. 
END  IF 
IF  (D.GT.3)  THEN 

FLAG(4)=.FALSE. 
END  IF 
IF(.N0T.(FLAG(1).AND.FLAG(2).AND.FLAG(3).AND.FLAG(4)))  THEN 

FLAG(1)=.TRUE. 

FLAG(2)=.TRUE. 

FLAG(3)=.TRUE. 

FLAG(4)=.TRUE. 

GO  TO   112 
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END  IF 

*  AFTER  SECOND  INSPECTION  DETERMINATION  OF  LCB 

LCB=SELCB(A,B,C,D) 

*  CHECKING  FOR  BONUS  AFTER  SECOND  INSPECTION 

IF  (  LCB.GT.BLINE  )  THEN 

BSI=BSI+1 

GO  TO   112 
END  IF 

*  COUNTING  FOR  CHECKING  1000  REPLICATIONS 
110     COUNT  =  COUNT  +1.0 

A=0.0 
B=0.0 
C=0.0 
D=0.0 

*  CHECKING  FOR  1000  REPLICATIONS 

IF  (COUNT. GT. 1000)  THEN 

GO  TO   130 
END  IF 
120   CONTINUE 

*  **AAAA-*rt*AA-AA-AAAAAAA-^^ 

*  PERCENTAGE  ESTIMATION  OF  GETING  BONUS  FOR  FIRM  A 
130   BTOT  =  BFI+BSI 

PRCT=BTOT/ 1000.0 

*  **************  A  A  A  A  A  A  A  A  *****  *****  A A**  A  A  A  AAA  A  A  A  A  A  A  ***************** 

*  WRITING  RESULTS 
WRITE  (15,1) 

1  FORMAT  (//,16X,*  BONUS  PLAN  SIMULATION  FOR  FIRM  A  * ,2X) 
WRITE  (15,2) 

2  FORMAT  (  1 6X ,  '  **  *  *  A  ♦'»  *  *  *  *  A  *  *  A  *  *  *  *  *  *  ************** '  2X ,  /  /  ) 
WRITE  (15,3) 
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3    F0RMAT(4X,'  FIRM  A  HAS  FOLLOWING  LONG  RUN  PROBABILITIES  IN  TESTS  ' 
+  ,2X) 
WRITE  (15,4) 

+./.2X) 
WRITE  (15,5)  PSIM,PSITH,PSIV,PSIA 

5  FORMAT(4X, 'PROBABILITY  OF  SUCCESS  IN  MANUFACTURER  TEST  IS',2X, 
+F8. 6, //,4X, 'PROBABILITY  OF  SUCCESS  IN  TEMP.  AND  HUM.  TEST  IS',2X, 
+F8. 6, //,4X, 'PROBABILITY  OF  SUCCESS  IN  VIBRATION  IS',2X, 

+F8. 6, //,4X, 'PROBABILITY  OF  SUCCESS  IN  ALTITUDE  IS ' ,F8 . 6,2X , //) 
WRITE  (15,6)  BFI 

6  F0RMAT(4X,*FIRM  A  GAT  BONUS  AFTER  FIRST  INSPECTION  ' ,2X,F6 . 1 ,2X , 
+'TIMES' ,2X) 

WRITE  (15,7) 

7  FORMAT  (  4X  ,  '  ^***********************************Vryr**yr*yr^**^yr*ycyr*^yr 

*******  \2X,/) 
WRITE  (15,8)  BSI 

8  F0RMAT(4X,*FIRM  A  GAT  BONUS  AFTER  SECOND  INSPECTION  ' ,2X,F6 . 1,2X, 
+'TIMES*,2X) 

WRITE  (15,9) 

9  FORMAT  (4X  '  *Mr***HHHMriHHHHr**f^^ 
+****** ' , 2X , / ) 

WRITE  (15,11)  BTOT 

11  FORM AT (4X,' TOT ALLY  FIRM  A  GAT  BONUS  IN  1000  REPLICATIONS  ',  2X, 
+F5.1,2X,*TIMES,,2X) 

WRITE  (15,12) 

12  FORM  AT  (4X  '  ***  *  *  ^H^*^^HM^^H^r*^Mr^W^r^^ 

+******', 2X,/) 

WRITE  (15,13)  PRCT 

13  F0RMAT(4X,'GETTING  BONUS  PERCENTAGE  OF  FIRM  A  IS' , IX, F6. 3 ,2X) 
WRITE  (15,14) 

14  FORMAT (4X,1 ' 
+******' ,2X,///) 

WRITE   (15,15)   BLINE 
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15  F0RMAT(4X,'  BONUS  LINE  FOR  FIRMS  IS ' , 1X,F6.4,2X) 
WRITE  (15,16) 

16  FORM AT (4X, '***** 
+*******, 2X,///) 

STOP 
END 
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APPENDIX  I.     95  %  LCB'S  FOR  DSBS  (  EQUAL  PROBABILITIES  ) 


Table   18.     95  °/J 

)  LCB'S  FOR  DOUBLE  SAMPLING  BONUS  SYSTEM 

FAILURE 
VECTOR 

95  %  LCB 

FAILURE 
VECTOR 

95  %  LCB 

(0  0  00) 

1.0000000 

(000  1) 

0.9099184 

(0  0  02) 

0.8706521 

(00  10) 

0.9375849 

(0011) 

0.8808635 

(00  12) 

0.8387057 

(0020) 

0.8968240 

(002  1) 

0.8528913 

(0022) 

0.8141473 

(0030) 

0.8717731 

(0031) 

0.8281843 

(0032) 

0.7881665 

(004  0) 

0.8422214 

(004  1) 

0.8021229 

(0  042) 

0.7645337 

(0  100) 

0.9099184 

(0101) 

0.8548708 

(0  102) 

0.8224728 

(01  10) 

0.8808635 

(0111) 

0.8368122 

(0  112) 

0.8010584 

(0  120) 

0.8548636 

(0121) 

0.8131180 

(0  122) 

0.7779465 

(0  130) 

0.8294836 

(0131) 

0.7867751 

(0  132) 

0.7541937 

(0  140) 

0.8001103 

(0141) 

0.7624319 

(0  142) 

0.7305543 

(02  00) 

0.8657608 

(020  1) 

0.8224728 

(0202) 

0.7905909 

(02  10) 

0.8387058 

(02  1  1  ) 

0.7996263 

(02  12) 

0.7653871 

(0220) 

0.8154084 

(022  1) 

0.7777317 

(0222) 

0.7414687 

(0230) 

0.7881665 

(0231) 

0.7541937 

(0232) 

0.7210189 

(0240) 

0.7624319 

(024  1) 

0.7326533 

(0242) 

0.6996729 

(1000) 

0.9945652 

(1001) 

0.9099184 

(1002) 

0.8706521 
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Table   19.     95  °A 

>  LCB'S  FOR  DOUBLE  SAMPLING  BONUS  SYSTEM 

FAILURE 
VECTOR 

95  %  LCB 

FAILURE 
VECTOR 

95  %  LCB 

(10  10) 

0.9375849 

(10  11) 

0.8808635 

(10  12) 

0.8369564 

(102  0) 

0.8968240 

(1021) 

0.8519021 

(102  2) 

0.8111921 

(10  30) 

0.8717731 

(1031) 

0.8258852 

(1032) 

0.7881665 

(104  0) 

0.8422214 

(1041) 

0.8001103 

(  1042) 

0.7631623 

(1100) 

0.9099184 

(110  1) 

0.8548708 

(1102) 

0.8224728 

(1110) 

0.8808635 

(1111) 

0.8368122 

(1112) 

0.8002886 

(1120) 

0.8528913 

(112  1) 

0.8121204 

(112  2) 

0.7768002 

(113  0) 

0.8258852 

(113  1) 

0.7862302 

(1132) 

0.7532073 

(1140) 

0.8001103 

(114  1) 

0.7616678 

(1142) 

0.7295684 

(  1200) 

0.8657608 

(12  0  1) 

0.8224728 

(1202) 

0.7872553 

(12  10) 

0.8387058 

(12  11) 

0.7987179 

(12  12) 

0.7626535 

(1220) 

0.8130434 

(1221) 

0.7758973 

(1222) 

0.7406143 

(12  30) 

0.7860773 

(1231) 

0.7538874 

(1232) 

0.7200024 

(124  0) 

0.7624319 

(1241) 

0.7309083 

(1242) 

0.6994067 

(2000) 

0.9891304 

(200  1) 

0.9099184 

(2002) 

0.8657608 

(20  10) 

0.9324048 

(2011) 

0.8808550 

(2  0  12) 

0.8341966 
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Table  20.     95  °A 

>  LCB'S  FOR  DOUBLE  SAMPLING  BONUS  SYSTEM 

FAILURE 
VECTOR 

95  %  LCB 

FAILURE 
VECTOR 

95  %  LCB 

(2020) 

0.8968240 

(2021) 

0.8499787 

(2022) 

0.8110244 

(2030) 

0.8717731 

(2031) 

0.8247706 

(2032) 

0.7845937 

(2  040) 

0.8422214 

(2041) 

0.8001103 

(2  042) 

0.7624320 

(2  100) 

0.9099184 

(2101) 

0.8548708 

(2  102) 

0.8224728 

(2110) 

0.8808550 

(2111) 

0.8362838 

(21 12) 

0.8003395 

(2120) 

0.8519021 

(2121) 

0.8121204 

(2  122) 

0.7750615 

(2130) 

0.8258852 

(2131) 

0.7846695 

(2  132) 

0.7517914 

(2  14  0) 

0.8001103 

(2141) 

0.7605843 

(2  142) 

0.7284731 

(22  00) 

0.8608695 

(220  1) 

0.8222825 

(22  02) 

0.7860733 

(22  10) 

0.8341965 

(22  11) 

0.7965909 

(22  12) 

0.7626535 

(222  0) 

0.8110244 

(222  1) 

0.7750614 

(2222) 

0.7401020 

(22  30) 

0.7845957 

(22  3  1) 

0.7519901 

(2232) 

0.7188347 

(2240) 

0.7611454 

(224  1) 

0.7295684 

(2242) 

0.6975686 
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APPENDIX  J.     DOUBLE  SAMPLING  BONUS  SYSTEM  WITH  EQUAL 

PROBABILITIES 

Table  21.     DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.825 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.825 

PROB/S 

BONUS  % 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.008 

0.9375 

0.030 

0.9450 

0.059 

0.9500 

0.079 

0.9600 

0.172 

0.9700 

0.351 

0.9750 

0.487 

0.9800 

0.632 

0.9850 

0.799 

0.9900 

0.920 

0.9950 

0.982 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  =  0.B25 
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Figure  5.      Double  Sampling  Bonus  System  With  LCBFB  =  0.825 
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Bonus  percentages  are  tabulated  and  plotted  below  with  LCBFB  0.850 


Table  22.     DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.850 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.850 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.007 

0.9375 

0.025 

0.9450 

0.050 

0.9500 

0.063 

0.9600 

0.144 

0.9700 

0.293 

0.9750 

0.423 

0.9800 

0.565 

0.9850 

0.738 

0.9900 

0.887 

0.9950 

0.970 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILmES;  LCBFB  =  0.B50 
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Figure  6.      Double  Sampling  Bonus  System  With  LCBFB  =  0.850 
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Bonus  percentages  are  tabulated  and  plotted  below  with  LCBFB  0.875 


Table  23.     DSBS  (EQUAL  PROBABILITIES)  LCBFB 

=  0.875 

LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.875 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.007 

0.9375 

0.023 

0.9450 

0.045 

0.9500 

0.057 

0.9600 

0.122 

0.9700 

0.236 

0.9750 

0.353 

0.9800 

0.487 

0.9850 

0.643 

0.9900 

0.799 

0.9950 

0.941 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  =  0.B75 
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Figure  7.      Double  Sampling  Bonus  System  With  LCBFB  =  0.875 
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Bonus  percentages  are  tabulated  and  plotted  below  with  LCBFB  0.900 


Table  24.     DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.900 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.900 

PROB/S 

BONUS  % 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.002 

0.9375 

0.010 

0.9450 

0.020 

0.9500 

0.027 

0.9600 

0.052 

0.9700 

0.138 

0.9750 

0.212 

0.9800 

0.319 

0.9850 

0.486 

0.9900 

0.671 

0.9950 

0.872 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  =  0.900 
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Figure  8.      Double  Sampling  Bonus  System  With  LCBFB  =  0.900 
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Bonus  percentages  are  tabulated  and  ploted  below  with  LCBFB  0.950 


Table  25.     DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.950 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.950 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.002 

0.9375 

0.010 

0.9450 

0.018 

0.9500 

0.022 

0.9600 

0.045 

0.9700 

0.107 

0.9750 

0.163 

0.9800 

0.237 

0.9850 

0.362 

0.9900 

0.502 

0.9950 

0.695 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILmES;  LCBFB  =  0.950 
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Figure  9.      Double  Sampling  Bonus  System  With  LCBFB  =  0.950 
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Bonus  percentages  are  tabulated  and  ploted  below  with  LCBFB  0.999 


Table  26.     DSBS  (EQUAL  PROBABILITIES)  LCBFB  =  0.999 


LOWER  CONFIDENCE  BOUND  FOR  BONUS  IS  0.999 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

PROB.'S 

BONUS  % 

0.9200 

0.002 

0.9375 

0.007 

0.9450 

0.010 

0.9500 

0.011 

0.9600 

0.028 

0.9700 

0.067 

0.9750 

0.106 

0.9800 

0.166 

0.9850 

0.274 

0.9900 

0.407 

0.9950 

0.626 

DOUBLE  SAMPLING  BONUS  SYSTEM 
EQUAL  TEST  PROBABILITIES;  LCBFB  =  0.999 
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Figure   10.      Double  Sampling  Bonus  System  With  LCBFB  =  0.999 


141 


APPENDIX  K.     95  %  LCB'S  FOR  DSBS  (DIFFERENT  PROBABILITIES) 

Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.825 

Table  27.     DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.825 


0.950 

0.975 

0.990 

0.995 

0.950 

0.079 

0.364 

0.632 

0.689 

0.975 

0.099 

0.487 

0.831 

0.882 

0.990 

0.109 

0.550 

0.920 

0.968 

0.995 

0.110 

0.558 

0.934 

0.982 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  =  0.825 
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Figure  11.      Double  Sampling  Bonus  System  With  LCBFB  =  0.825 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.850 

Table  28.     DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.850 


0.950 

0.975 

0.990 

0.995 

0.950 

0.063 

0.320 

0.612 

0.680 

0.975 

0.077 

0.423 

0.802 

0.871 

0.990 

0.083 

0.476 

0.887 

0.957 

0.995 

0.084 

0.483 

0.900 

0.970 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  =  0.850 
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Figure  12.      Double  Sampling  Bonus  System  With  LCBFB  =  0.850 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.875 

Table  29.     DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.875 


0.950 

0.975 

0.990 

0.995 

0.950 

0.057 

0.277 

0.566 

0.664 

0.975 

0.067 

0.353 

0.721 

0.846 

0.990 

0.072 

0.388 

0.799 

0.929 

0.995 

0.072 

0.392 

0.809 

0.941 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  =  0.875 


oSS  tests  p«ob.-s 


Figure   13.      Double  Sampling  Bonus  System  With  LCBFB  =  0.875 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.900 

Table  30.     DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.900 


0.950 

0.975 

0.990 

0.995 

0.950 

0.027 

0.164 

0.470 

0.605 

0.975 

0.031 

0.212 

0.603 

0.778 

0.990 

0.034 

0.235 

0.671 

0.860 

0.995 

0.034 

0.238 

0.679 

0.872 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  =  0.900 
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Figure  14.      Double  Sampling  Bonus  System  With  LCBFB  =  0.900 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.950 

Table  31.     DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.950 


0.950 

0.970 

0.990 

0.995 

0.950 

0.022 

0.126 

0.369 

0.506 

0.975 

0.029 

0.182 

0.502 

0.686 

0.990 

0.029 

0.182 

0.502 

0.686 

0.995 

0.029 

0.185 

0.509 

0.695 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  =  0.950 
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Figure  15.      Double  Sampling  Bonus  System  With  LCBFB  =  0.950 
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Bonus  percentages  are  tabulated  and  ploted  with  different  probabilities. 

LCBFB  =  0.999 

Table  32.     DSBS  (DIFFERENT  PROBABILITIES)  LCBFB  =  0.999 


0.950 

0.970 

0.990 

0.995 

0.950 

0.011 

0.057 

0.174 

0.229 

0.975 

0.021 

0.106 

0.302 

0.406 

0.990 

0.024 

0.151 

0.407 

0.570 

0.995 

0.026 

0.161 

0.463 

0.626 

DOUBLE  SAMPLING  BONUS  SYSTEM 
LCBFB  =  0.999 
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Figure  16.      Double  Sampling  Bonus  System  With  LCBFB  =  0.999 
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